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ABSTRACT 


Numerical  results  for  flow  and  heat  transfer  are  obtained  by 
finite  difference  method  for  steady  fully  developed  laminar  flow  in  curved 
rectangular  channels  under  the  thermal  boundary  conditions  of  axially 
uniform  wall  heat  flux  and  peripherally  uniform  wall  temperature  at  any 
axial  position. 

It  is  found  that  the  parameter  K  (Dean  number)  is  important  for 
the  governing  equations  under  consideration.  The  numerical  method  yields 
solutions  for  a  reasonably  wide  range  of  Dean  number  for  the  curved 
rectangular  channels  with  aspect  ratios  y  =  0.2,  0.5,  1,  2,  and  5  con¬ 
sidered  in  this  study.  The  method  used  complements  the  boundary  layer 
approximation  which  is  available  in  the  literature  and  is  applicable  for 
high  Dean  number  region  only. 

Graphical  results  for  Pr  =  0.73  are  presented  for  the  ratios 
of  Nusselt  numbers,  and  products  of  friction  factor  and  Reynolds  number 
between  the  curved  and  straight  channels,  Nu/(Nu)Q  and  (fRe)/(fRe)Q, 
respectively,  with  Dean  Number  as  a  parameter.  Typical  examples  for 
axial  velocity  distributions,  temperature  profiles,  isotherms,  and  stream¬ 
lines  and  velocity  profiles  for  secondary  flow  are  also  shown.  For  square 
channel,  the  effect  of  Prandtl  number  on  the  heat  transfer  results  is  also 
investigated. 

Comparison  of  the  present  result  with  the  results  for  the  curved 
square  channel  available  in  the  literature  shows  clearly  that  reasonable 
estimate  can  be  made  for  the  flow  and  heat  transfer  results  for  the  Dean 
number  ranging  from  150  to  1000  where  currently  accurate  solutions  are 
not  available. 
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NOMENCLATURE 


A 

a 
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C0,C1 ,C1 ' 
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G  r 
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h 

h 
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=  cross-section  area 
=  width  of  a  curved  rectangular  channel 
=  height  of  a  curved  rectangular  channel 
=  constants  defined  in  equation  (3-2) 

=  -  c-j  De3/4vy 

=  specific  heat  at  constant  pressure 

=  constant  defined  in  equation  (2-12) 

=  axial  pressure  gradient,  -  3PQ/3Z 

=  temperature  gradient,  3T/3Z  or  3T/R  3$ 

=  equivalent  hydraulic  diameter,  4A/S  =  2ab/(a+b) 
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=  Eckert  number,  Id  /c  (T-T,,) 

c  p  W 

=  function  defined  in  equation  (3-1) 
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=  friction  factor,  2  ^/(pW2) 

=  heat  generation 

.  2 

=  Grashof  number,  D  ^gBCg/v 

=  gravitational  acceleration 

=  dimensionless  grid  spacing,  a/MDg  =  b/2  NDg 

=  average  heat  transfer  coefficient 

1/2 

=  Dean  number,  Re(D  /R  )  '  ,  see  equation  (2-19)  or 

C 

(4-11) 

=  thermal  conductivity 
=  number  of  divisions  in  the  X-direction 
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s  Reynolds  number,  a  W/v 
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=  dimensionless  radius  of  curvature  of  a  curved  rectangular 
channel 

=  circumference  of  cross-section 
=  local  temperature 
=  mixed  mean  temperature  difference 
»  wall  temperature 
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in  cylindrical  coordinates 
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=  thermal  diffusivity 
=  coefficient  of  thermal  expansion 
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CHAPTER  I 


INTRODUCTION 

The  influence  of  the  dominating  body  force  on  the  flow  and 
heat  transfer  characteristics  in  pipes  or  channels  of  various  cross 
sectional  shapes  has  been  the  subject  of  many  investigations  in  recent 
years .  The  forced  convection  problems  with  body  force  effect  are  fre¬ 
quently  encountered  in  various  heat  exchangers,  cooling  or  heating 
systems,  reactors  and  heat  engines .  The  familiar  sources  of  body  force 
acting  on  fluid  element  are  gravitational,  centrifugal,  Coriolis,  mag¬ 
netic  and  any  artificial  acceleration  forces.  It  is  noted  that  body 
force  also  results  due  to  thermal  radiation  but  this  effect  is  negligibly 
small  in  most  problems  of  interest.  In  general  a  body  force  can  act  at 
any  angle  with  the  main  flow  in  the  channel.  However  in  this  study  con¬ 
sideration  will  be  given  only  to  the  case  where  the  body  force  is  acting 
normal  to  the  main  flow.  Due  to  the  body  forces  acting  in  the  cross 
section  perpendicular  to  the  main  flow,  the  flow  is  now  three-dimensional 
and  helical  even  for  the  fully  developed  situation.  This  means  that  the 
additional  two  velocity  components  arise  in  the  plane  of  the  cross  section 
and  these  velocity  components  together  are  usually  called  the  secondary 
flow. 

Under  certain  conditions,  a  secondary  flow  due  to  body  force 
will  not  be  established  until  a  critical  value  of  the  characteristic 
parameter  based  on  body  force  and  the  flow  field  is  reached.  The  stability 
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problem  relating  to  the  onset  of  the  secondary  flow  due  to  centrifugal 
or  Coriolis  forces  and  the  thermal -convective  instability  problem  (due 
to  buoyant  force)  have  also  attracted  much  attention  in  the  past  [1,2,3]. 

A  number  of  combinations  concerning  the  kinds  of  body  forces 
and  the  geometrical  shapes  of  the  flow  passages  are  possible  in  forced 
convective  heat  transfer  problems  [4,5],  However,  in  this  study  con¬ 
sideration  will  be  given  only  to  a  particular  combination  dealing  with 
centrifugal  force  and  rectangular  channel,  namely,  curved  rectangular 
channel.  For  a  flow  in  the  curved  channel,  the  critical  value  of  the 
characteristic  parameter  (Dean  number)  is  zero.  As  the  Dean  number  in¬ 
creases,  the  intensity  of  the  secondary  flow  increases;  consequently, 
the  flow  and  temperature  fields  are  affected  by  the  secondary  flow. 

Since  forced  convective  heat  transfer  in  curved  pipes  is  closely  related 
to  the  present  problem,  a  review  of  the  pertinent  works  on  flow  and  heat 
transfer  in  curved  pipes  is  in  order. 

1 . 1  Review  of  the  Literature  on  Flow  and  Heat  Transfer  in  Curved  Pipes 

The  literature  on  flow  in  curved  pipes  is  very  extensive  be¬ 
cause  of  its  technical  importance.  The  effect  of  centrifugal  force  on 
flow  in  pipe  bends  was  first  mentioned  by  Thomson  [6]  in  1876.  Experi¬ 
mental  investigations  of  the  secondary  flow  in  curved  pipes  were  made 
by  Eustice  [7,8],  White  [9],  Taylor  [10]  and  Yarnell  and  Nagler  [11]. 

The  theoretical  analysis  was  first  given  by  Dean  [12,13]  for  fully  de¬ 
veloped  laminar  flow  of  incompressible  fluid  in  a  curved  pipe  of  circular 
cross-section.  Dean  employed  a  perturbation  method  and  his  analysis  is 
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applicable  only  in  a  low  Dean  number  region  (0  <_  K  <  ~  25).  However, 
Dean's  work  appears  to  have  sparked  a  series  of  works  on  curved  pipes 
thereafter.  Adler  [14]  introduced  the  important  concept  of  boundary 
layer  for  secondary  flow  along  the  wall  in  his  analysis  for  laminar 
flow  with  large  Reynolds  number.  His  experimental  work  covers  a  wide 
range  encompassing  both  laminar  and  turbulent  flows  and  involves  the 
measurements  of  velocity  distributions.  Barua's  [15]  analysis  on 
secondary  flow  in  stationary  curved  pipes  employed  boundary  layer  approxi 
mation  along  the  wall  for  the  laminar  flow  when  the  Dean  number  is  large 
and  complements  Dean's  work.  Ito  [16]  proposed  accurate  empirical 
formulas  on  friction  factors  for  fully  developed  turbulent  flow  in  curved 
pipes  and  the  empirical  equation  for  transition  from  laminar  to  turbulent 
flow.  He  also  presented  asymptotic  expressions  for  friction  factors 
deduced  from  the  1/7  th-power  and  the  logarithmic  velocity  distribution 
laws  based  on  boundary  layer  approximation.  Ito  [17]  also  reported  theo¬ 
retical  and  experimental  investigation  on  laminar  flow  in  curved  pipes. 
Weske  [18]  found  that  the  turbulent  flow  in  curved  ducts  may  be  analyzed 
by  methods  adopted  from  the  theory  of  boundary  layers.  Cuming  [19]  gave 
solutions  for  flow  in  curved  pipes  of  circular,  elliptic  and  rectangular 
sections  in  power  series  of  the  curvature  of  the  pipe.  Topakoglu  [20] 
employed  similar  method  for  steady  laminar  flow  of  an  incompressible 
viscous  fluid  in  curved  pipes  of  circular  and  concentric  annular  sections 
By  assuming  that  the  actual  secondary  motion  is  approximately  a  uniform 
stream  for  the  most  part  of  the  section,  Dean  and  Hurst  [21]  obtained 
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solution  for  the  motion  of  fluid  in  curved  pipes  of  circular  and  square 
sections.  Hawthorne  [22]  and  Detra  [23]  analyzed  the  secondary  flow 
problems  in  curved  pipes  based  on  inviscid-fluid  theory.  It  is  seen 
that  the  published  works  on  flow  in  curved  pipes  are  numerous;  the 
above  review  deals  only  with  the  literature  which  is  readily  available. 

Jeschke's  experimental  work  on  heat  transfer  for  turbulent 
flow  of  air  in  helically  coiled  tubes  was  quoted  by  McAdams  [24]. 

Hawes'  [23]  experimental  investigation  concerned  with  the  measurements 
of  the  velocity  and  temperature  distributions  of  the  water  in  a  coiled 
pipe,  but  his  temperature  profiles  showed  somewhat  unusual  character. 
Pratt's  [26]  work  is  directed  toward  developing  empirical  equations  of 
turbulent  flow  for  both  the  internal  and  external  film  coefficients  of 
heat  transfer  in  a  cooling  coil  immersed  in  a  stirred  liquid  contained 
in  a  tank.  Ede  [27]  pointed  out  the  lack  of  heat  transfer  data  on 
various  pipe  bends  and  reported  the  results  on  local  heat  transfer  coef¬ 
ficient  for  water  in  a  right-angled  pipe  bend  covering  a  range  of  Reynolds 
number  from  500  to  50,000.  He  noted  that  with  laminar  flow  the  effect 
of  the  ratio  of  bend  radius  to  pipe  radius  is  large.  Recently,  Seban 
and  McLaughlin  [28]  presented  friction  and  heat  transfer  results  for  the 
laminar  flow  of  oil  and  the  turbulent  flow  of  water  in  tube  coils  with 
uniform  heating  having  the  ratios  of  coil  to  tube  diameter  of  17  and  104, 
for  Reynolds  numbers  from  12  to  65,000.  They  too  noted  the  need  for  ad¬ 
ditional  experimental  data  to  properly  define  the  heat  transfer  in  curved 
tubes.  Rogers  and  Mayhew  [29]  provided  additional  experimental  results 
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for  forced  convection  heat  transfer  and  friction  factors  using 
water  in  steam  heated  helically  coiled  tubes  with  turbulent  flow. 

At  this  point  it  is  quite  clear  that  additional  experimental 
investigations  and  theoretical  analysis  for  heat  transfer  and  pressure 
loss  in  curved  pipes  are  still  badly  needed.  Mori  and  Nakayama  [30,  31 , 
32]  presented  a  series  of  comprehensive  studies  on  forced  convective 
heat  transfer  in  curved  pipes  for  both  laminar  and  turbulent  regions. 

The  results  of  theoretical  analysis  based  on  boundary  layer  approximation 


along  the  pipe  wall  for  fully  developed  laminar  flow  in  a  curved  pipe 
under  the  condition  of  uniform  heat  flux  at  large  Dean  numbers  were 
shown  to  be  in  good  agreement  with  experimental  study  using  air.  For 


fully  developed  turbulent  flow  under  the  condition  of  uniform  heat  flux, 
the  results  of  theoretical  analysis  assuming  a  boundary  layer  along  the 
pipe  wall  and  experimental  results  using  air  were  again  shown  to  be  in 
good  agreement.  The  theoretical  analysis  under  the  condition  of  uniform 
wall  temperature  and  the  presentation  of  practical  formulae  for  both 
laminar  and  turbulent  flows  rounded  out  this  series  of  investigations. 
Kubair  and  Kuloor  [33]  presented  empirical  correlations  on  pressure  drop 

and  heat  transfer  for  laminar  flow  in  coiled  pipes  under  constant  tube 

•  • 

wall  temperature  conditions.  Very  recently,  Ozisik  and  Topakoglu  [34] 
presented  heat  transfer  results  for  hydrodynamical ly  and  thermally  fully 
developed  laminar  flow  in  a  curved  pipe  under  the  conditions  of  axially 
uniform  heat  flux  and  uniform  peripheral  wall  temperature  at  any  cross 


ft 


6 


section  for  small  Dean  number  region  using  a  method  of  series  expansion; 
the  results  are  applicable  for  small  curvatures.  It  is  noted  that  solu¬ 
tion  for  laminar  forced  convection  in  a  curved  pipe  for  the  intermediate 
Dean  number  region  is  still  unavailable  and  some  work  can  be  expected 
in  the  near  future  for  this  region.  Finally,  it  may  be  of  some  interest 
to  note  that  analytical  solutions  [35,36]  for  flow  of  some  non-Newtonian 
fluids  in  curved  pipes  are  also  available, 

1 o2  Review  of  the  Literature  on  Flow  and  Heat  Transfer  in  Curved 
Rectangular  or  Other  Noncircular  Channels 

In  contrast  to  the  number  of  published  papers  dealing  with 
flow  and  heat  transfer  characteristics  in  curved  pipes,  the  published 
works  relating  to  flow  and  heat  transfer  in  rectangular  or  other  non¬ 
circular  channels  are  rather  limited  in  spite  of  its  significance  in 
practical  appl ications . 

Ito  [37]  first  reported  theory  on  laminar  flows  through  curved 
pipes  of  elliptic  and  rectangular  cross  sections  using  series  expansion 
in  Dean  number.  As  noted  earlier  Cuming  [19]  also  presented  theoretical 
analysis  for  flow  in  curved  pipes  of  elliptic  and  rectangular  sections 
using  perturbation  method.  Dean  and  Hurst  [21]  as  mentioned  before  ob¬ 
tained  some  analytical  results  for  laminar  flow  in  curved  square  channel 
by  assuming  uniform  stream  for  secondary  flow.  Using  boundary  layer 
approximation  along  the  wall,  Ludwieg  [38]  presented  analytical  results 
for  friction  factors  for  fully  developed  laminar  flow  in  helically  coiled 
square  channels  rotating  around  its  axis.  Extensive  experimental  results 
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covering  both  laminar  and  turbulent  flows  are  compared  against  the  theo¬ 
retical  results  for  laminar  flow  only  with  good  agreement,  Experimental 
results  for  friction  factors  are  also  presented  for  the  case  of  station¬ 
ary  curved  square  channel  with  laminar  and  turbulent  flows,  Kapur, 

Tyagi  and  Srivastava  [39]  solved  the  fully  developed  laminar  flow  problem 
in  a  curved  annulus  of  concentric  circular  cross  sections  for  the  case 
when  the  radius  of  curvature  of  the  annulus  is  large  compared  with  the 
radius  of  the  outer  curved  pipe.  Topakoglu's  [20]  solution  for  secondary 
flow  between  two  concentric  torus  shaped  pipes  shows  the  existence  of 
four  vortices  as  compared  to  a  pair  of  vortices  found  by  Kapur  et  al. 
based  on  simplified  analysis,  Eichenberger  [40]  analysed  entrance  region 
problem  in  a  curved  rectangular  channel  with  secondary  flow  by  assuming 
an  inviscid  flow.  Experimental  investigations  on  fully  developed  turbulent 
flows  in  a  plane  curved  channel  between  concentric  circular  walls  were 
reported  by  Wattendorf  [41]  and  Eskinazi  and  Yeh  [42]. 

Using  boundary  layer  approximation,  Mori  and  Uchida  [43]  pre¬ 
sented  analytical  results  for  fully  developed  laminar  flow  in  a  curved 
square  channel  under  the  thermal  condition  of  axially  uniform  wall  temper¬ 
ature  gradient.  Their  results  for  flow  and  heat  transfer  are  applicable 
only  for  the  regime  where  Dean  number  is  large.  Velocity  and  temperature 
fields  were  obtained  by  dividing  the  cross  section  into  core  and  boundary 
regions  and  considering  the  balances  of  kinetic  energy  and  entropy  pro¬ 
duction  for  the  boundary  layers. 

Neglecting  secondary  flow  effect,  Kreith  [44]  calculated  the 
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Nusselt  numbers  for  turbulent  flow  in  a  curved  channel  formed  by  the 
annulus  between  two  concentric  cylinders  with  the  fluid  flowing  in  a 
direction  perpendicular  to  their  common  axis  showing  good  agreement  with 
the  results  from  heat-transfer  experiments.  Mori  and  Uchida  [45]  also 
carried  out  theoretical  investigation  on  forced  convective  heat  transfer 
for  fully  developed  laminar  and  turbulent  flows  in  a  curved  channel 
under  the  condition  of  constant  wall  heat  flux  and  compared  the  theo¬ 
retical  results  with  experimental  measurements.  Ustimenko,  Zhurgembaev 
and  Nusupbekova  [46]  presented  flow  and  heat  transfer  results  for  fully 
developed  laminar  flow  in  curved  flat  channels  with  different  ratios  of 
heat  fluxes  at  the  inner  and  outer  walls  with  no  secondary  flow  effect. 

It  is  quite  clear  from  the  above  brief  review  that  more  theo¬ 
retical  and  experimental  work  is  certainly  required  to  bring  our  knowledge 
of  flow  and  heat  transfer  in  curved  channels  comparable  to  the  level 
provided  by  the  rather  extensive  data  relating  to  flow  and  heat  transfer 
in  straight  channels. 

1 0 3  Purpose  of  Present  Investigation 

The  purpose  of  this  work  is  to  present  an  accurate  numerical 
solution  on  forced  convective  heat  transfer  for  steady  fully  developed 
laminar  flow  in  curved  rectangular  channels  with  various  aspect  ratios 
under  the  thermal  boundary  conditions  of  axially  uniform  wall  heat  flux 
and  peripherally  uniform  wall  temperature  at  any  axial  position.  It  is 
noted  that  Mori  and  Uchida's  [43]  analytical  solution  using  boundary  layer 
approximation  is  valid  only  for  the  flow  regime  when  the  Dean  number  is 


■ 


„  •  '  .  • 


9 


large.  The  numerical  method  complements  the  boundary  layer  technique  for 
high  Dean  number  regime  and  yields  accurate  solutions  up  to  a  reasonably 
high  Dean  number  region  for  the  aspect  ratios  y  -  0,2,  0.5,  1,  2  and  5 
considered  in  this  study.  This  work  is  confined  to  a  study  of  the  effect 
of  centrifugal  force  field  on  fully  developed  laminar  flow  and  heat 
transfer  characteristi cs  in  curved  rectangular  channels  assuming  constant 
physical  properties;  it  deals  only  with  one  type  of  body  force  and  one 
kind  of  geometrical  shape.  Moreover  this  study  concerns  only  with  one 
aspect  of  a  broad  class  of  forced  convection  heat  transfer  problems  in 
rotating,  bent  and  heated  horizontal  tubes  involving  secondary  flow  caused 
by  body  forces.  In  order  to  understand  the  physical  meaning  of  the  various 
dimensionless  parameters  pertaining  to  the  problem  under  consideration 
and  the  relative  importance  of  each  term  in  the  differential  equations 
expressing  the  conservation  of  momentum  and  energy  for  an  incompressible 
fluid  in  a  curved  channel,  a  general  analysis  for  the  continuity,  Navier- 
Stokes  and  energy  equations  will  be  made  next. 


CHAPTER  II 


THEORETICAL  ANALYSIS 

2 . 1  General  Analysis  of  Fundamental  Equations  for  the  Curved  Channel 

Flow 

Consider  a  laminar  flow  with  constant  viscosity  and  thermal 
conductivity  coefficients  in  a  channel  of  rectangular  cross-section  whose 
horizontal  axis  is  bent  into  a  circle  of  radius  R  with  buoyancy  force 
effect.  Taking  the  origin  of  cylindrical  coordinates  (R,  $,  Y)  at  the 
center  of  curvature  with  (R,  $)  in  the  horizontal  plane,  the  differential 
equations  expressing  the  conservation  of  mass,  momentum  and  energy  for 
a  viscous  heat  generating  fluid  are 
Continuity  Equation: 


3U  .  U  1  3W  .  3V 
3R  T  R  R  3$  1’  3Y 


(2-1) 


R-Momentum  Equation: 
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^-Momentum  Equation: 
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Y-Momentum  Equation: 
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where  P  =  P'(R,Y)  +  Pq($)  -  3gY 


Energy  Equation: 


U9T  W9X  V9T 
U  9R  R  9$  V  9Y 
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(2-5) 


Introducing  representative  lengths,  a  and  b,  of  the  rectangular 
cross-section  of  the  channel,  representative  radius  of  curvature,  R  ,  of 
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the  channel,  representative  angle  ft  for  $,  and  the  reference  velocity 
components  U  ,  V  and  W  in  the  R,  Y  and  $-di recti ons ,  respectively, 
the  following  non-dimensional  quantities  may  be  defined: 


R  =  R  +  ax  =  R  ( 1  +  Ax)  ,  X  =  ax  , 

C  C  cl 


Y  =  by,$  =  ftct>,U  =  Uu,V  =  V  v,  and  W  =  W  w  . 


Using  the  above  relations,  the  continuity  equation  (2-1)  is 
reduced  to  the  following  non-dimensional  form. 


U 


/_C\  (h  J-2MJ.  ,  ,  u 

lV  }  laj  i9x  Aa  (1  +  Aax) 


,  0-K  An  1  .  9v_  =  n 

Q  a  U  (1  +  A  x)  94>}  9y  "  0 
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(2-6) 


since  9u/9x  and  9v/9y  may  be  expected  to  be  of  the  same  order  of  magni 
tude,  one  may  set 


ur  h 

(/)  (f)  -  0[1] 


(2-7) 


The  continuity  equation  now  becomes 
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Eliminating  pressure  P'  from  equations  (2-2)  and  (2-4)  and 
introducing  the  vorticity  function  c  =  9v/9x  -  9u/9y,  the  following 
relation  can  be  obtained  by  considering  that  the  centrifugal  force  term 
and  the  inertia  terms  are  equally  significant 


U 


c 


1 


(2-9) 


This  signifies  the  importance  of  the  secondary  motion.  The  resulting 
momentum  equation  in  non-dimensional  form  for  the  secondary  flow  is 


-1  9 


1  1 
-1  "2  ,  7  1 


Y  9x  +  9y  y  Aa  ^  +  x  x) 

a 

1  .1 

/  -1  9w  9v  9w  9U\  ,  ^-1  yZ  Z  w  9£ 

(y  37  3?  -  Y  37  3?>  ^  XY  (1  *  Ax)  3^ 

a 


i  i 


D  -1  2  .  2  rn  2r  .  2  1 

=  Rea  y  Xa  [Vy  c  +  a 


_  a2c 

a  (1  +  Ax)  2 


a  7  9<J> 


1  1  2 

,  r\  2  ^  2  1  9  w  ,  ,  1  9C 

^  aAa  (1  +  Ax)  94>9y  Aa  (1  +  \x)  9x 


^  ,  2  1  9 u -i  ,  r  -  -2,  -1  -1  90  1 

Aa  (1  +  Ax)2  9y  GrRea  Aa  Y  9x  "  Y  (1  +  Ax) 


9(w) 

9y 


(2-10) 


2  2 
„  2  9,29 

VY  =  ~  +  Y  “1 

Y  9x^  9y 


where 
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cr  = 


a  ftR 


Similarly,  the  axial  momentum  equation  in  non-dimensional  form  becomes 


_1  1  _1 

vi,  ±  3w  x  p,  2  2  2 r  w  3w-| 

Yu  37  +  v  37  +  n  °a  Y  C-rr  TTY) 

a 


1  1 


_  D  “1  •>  2  2  r  a _ 1  3Pin  2 

R  a  Aa  Y  yW  Q  (l  +  X  x)  3<f>  Vy  W 

C  a 


+  a 


2  1  1 

1  3  w  ,  2,  2  2 

7  7T  +  2oa  Aa  Y 


(1  +  X  x)  3(f)* 


1  3_u 

(1  +  X  x)  3<f) 

a 


,  1  3w  ^  2  W  -|  uw 

a  (1  +  Aax)  3x  -  Aa  (1  +  x  x)2^  '  Aa  HTx^I 


(2-11) 


The  reference  velocity  W  may  be  regarded  as  the  average  velocity  in  a 

V 

straight  channel  given  as 


W 


c 


_a_  3P  J_ 

y^  3(j>  Cq 


(2-12) 


where  Cq  is  a  constant  depending  on 
section. 

For  the  forced  convection 
gradient  is  usually  known  or  given, 
ature  difference  may  be  defined  as 


the  geometrical  shape  of  the  cross¬ 
problem,  the  axial  temperature 
Consequently,  the  reference  temper- 


T  -  T.,  =  ©  0 

w  c 
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It  is  clear  that  the  convective  term  involving  the  axial  temper¬ 
ature  gradient  and  one  lateral  conduction  terms  in  the  energy  equation 
must  be  equally  significant.  From  this  one  may  define  a  characteristic 
temperature  difference  as  follows: 


0 


c 


W  c  a 
c  2 

a 


(2-13) 


Then,  the  non-dimensional ized  energy  equation  neglecting  heat  generation 
is  given  as 


99  ■  w  99  ■  99 

9x  (1  +  X  x)  9c|)  v  9y 


1  1 


1  1 


Rea_1  Pr"1  [Xa2X  2V^2e  +  X  2  1 


(1  +  Xax) 


90_ 

9x 


1  9  0-]  r ^  n ^  ”1  ro„  r/9U\2  /9v^2-| 

+  °a  (1  +  Xax')"^7]  +  EC  Rea  {2y[(37}  +  1 


+  2A  - 2 - 2  (n_1  fjf  +  Xa2  u)2 

a  (1  +  Ax)2  3<t>  a 

a 


1 1 

,  ,  r,  -1.  3w  .  ,  '2  '2 
+  Xa[Xa  X  W  Xa  Y 


1  9w-i2-> 

(1  +  X  x)  9(j)J  i 


ht»  sriJ  gntvrovnt  msi  avfioavnoo  srti  69 

.  _  r— *..r  ^  hn«  trt<ir 
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S  >$  ,"5  :  \q  1 
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1  1 


-ii  9Un2 


,  r,  2  "2  9w  ,  "2  2  |£] 

+  Xa  [Xa  Y  “  +  X*  Y  9Y 


9x 


_1  1 

_i_  ^  r *\  2  2  ^”1  1  9u  ,  ^  -I  9w  W  n 2 

VAa  Y  a  (1  +  A  x)  3<t>  xa  3x  '  (1  +  A  x)^  }’ 


(2-14) 


For  a  fully  developed  flow  and  neglecting  viscous  dissipation,  one  may 
write 


Ec  «<  0  [1] 

n  >»  o  [i]  . 

The  last  condition  implies 

o  =  A  /ft  <«  0  [1] 

d  a 

since  geometrically  Aa  is  restricted  as 

A  <<  2  . 

a 

Define  the  following  quantities 

A  =  De/Rc  ,  a  =  De/Rcfl  «<  0  [1] 


and  y  =  a/b  =  0  [1]  where  Dg  =  2  ab/(a  +  b) 


I  I  's  L 
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and  replacing  Aa  and  a  in  equations  (2-8),  (2-10),  (2-11)  and  (2-14), 

a  a 

by  A  and  a  ,  respectively,  the  following  simplified  equations  result. 
The  Continuity  Equation: 


M  +  11=  0 

3x  3y 


(2-15) 


The  Momentum  Equation  for  Secondary  Flow: 


u^  +  v|t 

3x  3y 


1  1 

_1  “?  9  0 

Re  [A  VC  +  A 


3 

2 


1  3u 


t1  +  Ax)  +x  (i  +  xx)2 


] 


+  Gr  Re'2  A-1  -  1 


3(w)' 

3x  "  (1  +  Ax)  3y 


(2-16) 


The  Axial  Momentum  Equation: 


-1  c 

3w  ,  . .  3 w  _  n  ~ ~  1  r-,  2  0 

U  37  +  V  37  -  Re  (1  +  XX) 


.1  1 

.  , ”2  „2  ,  , 2  1  3w 
+  X  V  w  +  A  {1  +  Xxj  37 


-  A 


w 


(1  +  Ax)' 


r]  -  A 


uw 


(1  +  Ax) 


(2-17) 


The  Energy  Equation: 


36  ,  90  ,  w  30_ 

u  3x  v  3y  (1  +  Ax)  3<j) 


1 


1 


D  -1D  -1  r  "2  n2o  x  \ 2  1  30-i 

Re  Pr  [A  V  6  A  +  x"x )"  3x^ 


(2-18) 
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where  V2  = 

9x  9y 

It  may  be  of  interest  to  note  that  the  last  term  on  the  right- 
hand  side  of  equation  (2-16)  represents  a  centrifugal  force  effect  due 

to  curvature  and  the  term  next  to  the  last  represents  buoyancy  force 

*?• 

effect.  It  is  also  seen  that  the  last  term  on  the  right-hand  side  of 

equation  (2-17)  represents  Coriolis  force  effect. 

As  stated  earlier,  for  a  fully  developed  flow  condition 

Q  »>  0[1],  and  this  also  implies  that  X  «  0[1].  Thus  all  the  terms 
1/2  1  3/2 

involving  A  ,  A  and  X  may  be  neglected  as  compared  with  the  term 

-1/2 

of  order  one.  Similarly,  the  term  involving  X  must  be  retained  by 

considering  the  contribution  due  to  other  non-dimensional  characteristic 

parameters.  It  is  noted  that  Topakoglu  [20]  and  Ozisik  and  Topakoglu  [34] 

considered  the  cases  of  small  curvatures  for  laminar  flow  and  heat  transfer 

in  a  curved  pipe;  however,  the  results  of  analysis  show  that  the  terms 
1/2  3/2 

involving  A  ,  A  and  X  seem  to  have  negligible  effect  on  the  flow  and 
heat  transfer  results  at  least  for  the  very  low  Dean  number  region. 

The  following  general  conclusions  may  be  drawn  from  the  order 
of  magnitude  analysis: 

(1)  Re  =  0[lCf2]  and  X  =  0[UT2] 

Observations  of  the  equations  of  motion  show  that  inertia  terms, 
centrifugal  force  term  and  Coriolis  force  term  can  be  neglected  for  this 
case.  In  addition,  the  buoyancy  force  term  can  be  neglected  only  if 

r  _  A 

Gr  <  0[10  ]  and  must  be  considered  if  Gr  >  0 [1 0  ].  Observation  of  the 

energy  equation  shows  that  convective  terms  may  be  neglected  if  Pr  <  0 [1 0] . 


■ 
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3 

On  the  other  hand,  convective  terms  cannot  be  neglected  when  Pr  >  0 [1 0  ] . 
(2)  Re  >  0[10] 

In  contrast  to  the  above  case,  the  inertia  terms  and  the  con¬ 
vective  terms  must  be  considered  for  this  case.  The  buoyancy  force 

1  /2 

term  now  depends  on  the  order  of  magnitude  for  Re  A  .  The  following 
case  may  be  of  interest: 

1 

Re  A2  >  0 [ 1 02]  . 

3 

It  is  noted  that  the  buoyancy  force  term  must  be  considered  if  Gr  >  0 [1 0  ]. 
It  is  very  likely  that  this  case  may  occur  in  practice.  On  the  other 
hand,  if  Gr  <  0 [1 0  ],  the  buoyancy  force  term  can  be  neglected.  This 
condition  corresponds  to  the  problem  under  study.  It  is  also  noted 
that  the  Coriolis  force  effect  is  important  if  A  =  oCl].  However,  this 
condition  hardly  ever  occurs  for  the  fully  developed  flow  situation. 

The  foregoing  general  analysis  shows  clearly  that  the  following 
dimensionless  parameters  are  of  importance  for  flow  and  heat  transfer 
characteristics  in  a  curved  channel: 


Dimensionless  parameter  representing 

curvature  effect 

,  De 
"  Rc 

Reynolds  number 

Wc 

Re  =  D  — 
e  v 

Prandtl  number 

Pr  =  — 

a 

4  C2 
Gr  =  DeV  -f 


Grashof  number 


.  por]o  <  < 


0 


t 


-  -  Aa 


. 


■  ■ 


Observation  of  the  momentum  and  energy  equations  shows  that 


the  Reynolds  number.  Re,  and  the  curvature  effect  parameter,  A,  appear 

1  /2 

together  as  Re  A  =  K.  This  last  characteristic  parameter  for  flow 

in  a  curved  channel  is  usually  referred  to  as  Dean  number  in  the  liter¬ 
ature. 


In  order  to  gain  physical  insight  into  the  meaning  of  this 
important  characteristic  parameter,  the  Dean  number,  may  be  rewritten 
as  [2] 

2  2 
Wc  Wc 

pd~  prT~  1/2  1/2 

K  ■  [  if-  •  -T7 —  ]  =  (Re  •  r)  (2-19) 

c  c 


One  sees  that  the  Dean  number  consists  of  the  product  of  two  non-dimensional 
characteristic  parameters  Re  and  r.  The  meaning  of  Reynolds  number  is 
well  understood.  The  parameter  r  represents  the  ratio  of  the  centrifugal 
force  effect  to  the  viscous  force  effect  and  might  be  called  as  the 
centrifugal  Reynolds  number.  Finally,  one  notes  that  the  heat  generating 
term  in  the  energy  equation  (2-5)  was  included  to  show  that  the  analysis 
can  be  carried  out  for  heat  generating  fluid.  However,  for  simplicity, 
the  heat  generation  will  be  omitted  in  the  remaining  analysis. 

2.2  Basic  Equations  for  Fully  Developed  Flows  in  Curved  Rectangular 

Channels  and  Boundary  Conditions 


Consider  a  steady  fully  developed  laminar  flow  of  viscous  fluid 
in  a  curved  rectangular  channel  under  the  thermal  boundary  conditions  of 
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axially  uniform  wall  heat  flux  and  peripherally  uniform  wall  temperature 
at  any  axial  position.  Because  of  the  non-linearity  of  the  basic 
equations  governing  the  problem,  the  following  simplifying  assumptions 
are  made  to  facilitate  the  analysis: 

1.  Velocity  and  temperature  fields  are  fully  developed. 

2.  The  only  body  force  is  a  centrifugal  force  caused  by  the 
curvature  of  the  channel. 

3.  The  radius  of  curvature  of  the  rectangular  channel  is  large 
compared  with  the  hydraulic  diameter  of  the  cross-section 
of  the  channel . 

4.  No  internal  heat  generation. 

5.  Physical  properties  are  constant. 

6.  Viscous  dissipation  is  negligible. 

Taking  the  origin  of  the  rectangular  coordinates  (X,Y,Z)  at 
the  center  of  the  rectangular  cross-section  as  shown  in  Figure  1  and 
applying  the  assumptions  stated  above,  the  governing  equations  for  the 
present  problem  are 

Continuity  Equation: 


(2-20) 


Momentum  Equation  for  Secondary  Flow: 


u  § +  v  If  ■  ^  - 


2 


c 


(2-21) 
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Axial  Momentum  Equation: 


ii  9W  +  v  ^ 
U  3X  +  V  W 


1^0  2 

p  IT  +  W 


(2-22) 


Energy  Equation: 


ii  ST  \/  3T  .  ,,  3T  _  n2T 
U  3X  +  V  3Y  +  W  3Z  aV  T 


(2-23) 


where 


P  =  8V;  9U_ 

^  8X  “  9Y 


(2-24) 


and  the  pressure  at  any  point  consists  of  two  parts. 


P  =  P0(Z)  +  P  1  (X  ,Y) 


(2-25) 


Boundary  Conditions: 


U=V=W=T-Tw=0 


along  upper  and  lower  walls, 
Y  -  ^  or  Y  -  -  ^  and  -  j 
<  X  <  | 


V 


9U  _  9T 
9Y  9Y 


0 


along  center  line,  Y  =  0  (2-26) 


and  -  2—^—2 


U  =  V  =  W  =  T-Tw  =  0 


along  inner  and  outer  walls, 

X  =  -|-orX=|-  and  -  j  £  Y  £ 


ro|  cr 


,,  D  I 


—  i 
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Introducing  the  following  transformations, 

X  =  Dex  ,  Rc  =  Derc  ,  Y  =  Dey  ,  U  =  (v/D0)  u  , 

V  =  (v/De)v  ,  W  =  (vc/De)w  ,  c  =  -  c1De3/4vy, 

T  -  Tw  =  (c2DePrc)6 


where 


2ab 

(a  +  b) 


3T 

3Z  • 


and  a  dimensionless  stream  function  ip 


u  = 


v  = 


Ml 

9y 


_ 

3x 


(2-27) 


the  above  equations  may  be  restated  in  the  following  dimensionless  forms. 
Momentum  Equation  for  Secondary  Flow: 


3ip  3 
3y  3x 


3 ±  _3_ 
3x  3y 


(2-28) 
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Axial  Momentum  Equation: 


9t*  9w  dip  9w 
dy  dx  ”  9x  dy 


V2w 


+  4 


(2-29) 


Energy  Equation: 


(9£96  _  9ia0}  =  y20  . 
v9y  9x  9x  9y 


(2-30) 


Because  of  symmetry  with  respect  to  the  X-axis,  it  is  only 
required  to  consider  the  lower  half  of  the  rectangular  cross-section 
(see  Figure  1).  Consequently,  the  boundary  conditions  are: 


along  lower  wal  1 ,  y  =  - 


2D. 


and  - 


<  x  < 


2D  -  -  2D 

e  e 


dip  _  9  ip  96  _  q 

ax  -  ay2  -  3y  - 


along  center  line,  y  =  0 


and  .  _a_<  x  <  ^ 

e  e 


(2-31) 


|i  =  |i=  w  =  e  =  o 

3x  9y 


along  inner  and  outer  walls 

and 


a  a 

x  =  -  —  or  x  = 


2D 


2D. 


.1  y  i  o 


2D 


The  above  set  of  equations  constitutes  a  formal  mathematical 


statement  of  the  problem  under  consideration.  In  contrast  to  the  forced 
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convection  problem  with  buoyancy  effect,  one  notes  that  for  the  problem 
at  hand,  the  momentum  equations  and  the  energy  equation  are  uncoupled; 
namely,  velocity  field  can  be  found  independently.  It  is  possible  to 
solve  the  above  set  of  equations  analytically  by  a  method  of  successive 
approximation  in  principle  as  demonstrated  in  the  literature  for  similar 
problems  [19,20,37],  but  the  process  is  very  tedious  and  the  solution 
quickly  diverges  with  the  increase  of  the  Dean  number.  In  view  of  the 
considerable  difficulties  with  the  analytical  method,  the  numerical 
solutions  of  the  exact  equations  using  converging  iterative  procedures 
will  be  presented  in  this  study.  The  numerical  solutions  for  a  similar 
set  of  differential  equations  having  non-linear  terms  are  available  in 
the  literature  [47  to  53]. 
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CHAPTER  III 


NUMERICAL  SOLUTION 

3.1  Finite-Difference  Approximation 

In  recent  years,  finite-difference  methods  have  been  shown 
to  be  a  powerful  tool  for  the  solution  of  natural  convection  problems. 
The  prediction  of  natural  convection  by  numerical  method  was  recently 
reviewed  by  Churchill  [47].  Wilkes  and  Churchill  [48]  made  a  study  of 
natural  convection  in  a  rectangular  enclosure  with  one  vertical  wall 
heated  and  the  other  cooled  by  an  implicit  alternating  direction  finite- 
difference  method.  Finite-difference  methods  were  also  used  by  Samuels 
and  Churchill  [49]  to  compute  hydrodynamic  instability  due  to  natural 
convection  in  an  enclosed  horizontal  rectangular  region  heated  from 
below.  Barakat  and  Clark  [50]  made  a  study  of  the  two-dimensional 
transient  laminar  natural  convection  in  a  partially  filled  liquid  con¬ 
tainer  using  an  explicit  finite-difference  technique.  Recently,  Hwang 
[51,52]  presented  numerical  solution  for  fully  developed  combined  free 
and  forced  laminar  convection  in  horizontal  rectangular  channels  using 
point  successive-overrelaxation  method.  DeVahl  Davis  [53]  also  reported 
numerical  solution  for  laminar  natural  convection  in  an  enclosed  rec¬ 
tangular  cavity  by  similar  success! ve-overrelaxation  method.  In  this 
study,  the  point  successive-overrelaxation  method  was  employed  to  solve 
a  set  of  elliptic  partial  differential  equations  and  the  associated 
boundary  conditions. 
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Consider  equation  (2-29)  which  may  be  regarded  as 


„2  _  r  /.  3wn 

V  w  F-j  (u  ,v ,  9x,  9y) 


(3-1) 


The  Laplacian  function  V  w  at  a  mesh  point  (i,j)  may  be  expressed  by 


7  wi,j  =  C0  wi,j  +  Cl(wi+l,j  +  +  Cl(wi,j+1  +  wi,j-l)  + 


(3-2) 


Using  Taylor  series  expansion,  the  following  equations  can  be  written 
at  the  four  neighboring  points  of  the  mesh  point  (i,j). 


,  (1)  .  h2  (2)  .  h3  (3) 

wi+l  ,j  wi , j  hwi,jx  2!  wi  ,jx  3!  wi  ,jx 

+  hlw  (4)  +  Hi  w  (5)  +  h!w  (6)  + 

4!  wi  ,jx  5!  wi  ,jx  6!  wi,jx  •" 

.  (1)  .  h2  (2)  h3  (3) 

W-  n  .  =  w.  .  -  hw.  .  +  7T7  W.  .  -  OT  w.  . 

i-l  ,j  i  ,j  i  ,JX  2!  1,JX  3T  1,JX 


+  £w.  <4)  -£w.  (.5)  +  £w.  <6>  -  . 
4!  i,jx  5!  i,jx  61  i  ,jx 

,  (1)  .  h2  (2)  .  h3 

W.  .  .  ,  =  w.  .  +  hw.  ...  +  o!-W.  /  +  TT  w  . 


(3) 

’i .  j+1  "i  >j  '  ""i,jy  2!  "i  ,jy  3!  "i  ,jy 

+  hlw  (4)  +  hiw  (5)  +  h£  w  (6) 

+  4!  wi  ,jy  5!  wi , jy  6!  wi,jy  +  •" 

.  (1)  ^  h2  (2)  h3 

wi,j-l  =  wi,j  -  hwi,jy  +  2Twi,jy  -  3Tw}?iy 

+  £w.  (4)  .  hiw.  (5)  +  hjw.  (6)  . 

4!  l  ,jy  5!  l  ,jy  6!  l  ,jy 


(3-3) 


■ 


' 

I 
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It  is  noted  that  identical  dimensionless  mesh  size  (h)  is  used  in  X  and 
Y-directions .  Substituting  equation  (3-3)  into  equation  (3-2)  and  evalu¬ 
ating  the  coefficients  Cq,  C-| ,  C-j ' ,  the  Laplacian  function  at  the  mesh 
point  (i,j)  becomes 


2  4  1 

V  w.  •  = - 9  w.  .  +  -5-(w.,,  .  +  w.  .  +  w.  +  w.  .  -,) 

1  *J  h2  1  »J  h2v  1+1  ,3  l-l, J  1  ,J+1  l  »j-l 


h2,  (4)  ^  (4h  , 

+  I2(wi ,jx  +  wi ,jy  )  +  ' ’ ' 


(3-4) 


Comparing  equations  (2-29)  and  (3-1),  the  function  F  at  the  mesh  point 
(i,j)  may  be  expressed  as 


Fi,j  =  2h(wi+l,j  '  wi-l  ,j ^  ui,j  +  2h(wi,j+l  ‘  vi  ,j 


4  '  6  (wi  jx*  ui,j  +  wi  Jy5  vi,j)  +  ••• 


(3-5) 


Equating  equations  (3-4)  and  (3-5)  and  neglecting  higher  order  terms, 
the  finite-difference  equation  for  the  axial  momentum  equation  (2-29) 
may  be  expressed  as 


wi,j  =  I(Vl,j  +  wi-l,J  +  wi  ,j+l  +  wi  ,j-l} 


h2  ,  wi+l  „i  ~  wi-1,j  |  ..  wi  ,j+1  '  wi,j-l 

2h  vi,j  2h 


(3-6) 


4} 


b„6  X  nt  b92U  *  (H)  est*  ritM  fflofansmib  r6Drinebr  tett  beJan  >t  I 


't-t,fw  f+t.'  t.t-r  t.  +t  ^  t.f 

. 


29i'oo9d  (t.f  i  ifoq 


' 


'  H  &|  '  \ '  wU  -  maw  1 1  m 


Applying  similar  method,  the  finite-difference  equation  for  the  energy 
equation  (2-30)  becomes 
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9.  .  =  i(  0 .  ,  •  +  0  •  i  .  +  0.  . ,  +  0 .  *i) 
i  ,j  4V  i+i  ,j  i-i  ,j  i  ,j+i  i  ,j-r 


0 .  .  -  0 . 

i+i..]  i-i ,.i  + 


(3-7) 


6  "6 

-  Tf  {PrCUi  i  ’'Sh  +  ^  i  "■1>,j  +  1Oh  1-1  ,J]  +  W.  .} 

4  u  1  ,j  2h  i  ,j  2h  1  ,J 


To  obtain  the  finite-difference  equation  for  equation  (2-28), 
one  notes  that  equation  (2-28)  may  be  regarded  as  the  inhomogeneous 
bi harmonic  equation  in  the  following  form. 


-  F2[(lx*|y)**  (fx*fy)vZ<,,W*  <3'8) 

Omitting  details,  the  finite-difference  equation  for  the  secondary  flow 
stream  function  \p  at  the  mesh  point  (i,j)  may  be  obtained  by  double 
application  of  the  procedure  transforming  equation  (3-1)  into  its  finite- 
difference  form.  The  result  is 


f-i.j  =  (H)  ‘  10  ui,j)  +1+1  .j  +  (lo  +  10  ui,j)  +i-l  ,j 


+  (777  -  £77  V.  .)  Ip.  •  .  1  +  (777  +  777  V.  .)  l|j.  . 

MO  10  i,j'  M  ,J+1  MO  10  l,j'  M.J-1 


+  20+  40  ui,j*  +i+2 +  20  '  40  ui.j)  +1-2 ,j 


!  ’*  ' 


r  ,0 


zuosnseomortnt  sill  26  b^gei  ^  tm  (8S'S>  nctisups  w-  -w  s"° 
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+  20  +  40  vi  ,j *  ’h  ,j+2  +  20  '  40  vi,j*  ,j-2 


+  10  +  40  ui,j  +  40  vi  ,j*  ’N+l.j+l 


(3-9) 


+  10  +  40  ui,j  "  40  vi,j*  ’h+l.j-l 


+  10  '  40  ui,j  +  40  vi,j*  Vl,j+1 


1  h  h 

+  TO  ‘  40  ui  ,j  ‘  40  vi  ,j*  ’N-l  ,j-l 


h3  c2 

120  ^  wi,j  (-  wi  ,j+2  +  8  wi,j+l  -  8  wi,j-l  +  »U.2) 


It  is  noted  that  for  the  derivative  3w/9y  in  equation  (2-28),  five-point 
formula  is  used. 

It  is  well  to  note  that  the  finite-difference  expression  for 
the  stream  function  at  the  mesh  point  next  to  the  boundary  takes  the 
following  forms  after  satisfying  the  boundary  condition  for  ip: 

At  the  points  i  =  2,  j  =  3,4  ...  N  +  1 ,  .  is 

^  sj 


^.j  ‘  ^3,j  (?0  ■  u2,j  10*  +  ^2 .j+1  (!o  ■  v2,j  10* 


i  ( 4  h  \ 

^ .0-1  %  v2,j  To* 


1  h  h 

+  ^.j+l  10  +  u2,j  40  +  v2,j  40* 


f+t,r+f*  (t,rv  r  +  Uu  w*  or -)  + 


t.sv+us  e,s“nT  ->  rn-e*  + 
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/  1  h  h  \ 

v3 ,j-l  '■  10  u2,j  40  "  v2,j  40 


+  1  /  1  h  \ 

V2,j-2  ''  20  ‘  v2,j  40' 


+  ^.j  ^  20  +  u2,j  40^ 


+  ^2  ,j+2  20  +  v2,j  40^ 


h3  c2 

T20  <77>  w2,j  {'  W2,j+2  +  8w2,j+l  '  *2,j-l  +  w2,j-2)]/ 


r  -i  1  h  -I 

u  "  70  '  u2,j  40J 


(3-10) 


At  the  points  i  =  M,  j  =  3,4  ...  N  +  1,  ^  .  is 


^M,j  tyl-l.j  ^10  +  UM, j  10^  +  ^M.j+l  *10  '  vM,j  10^ 


+  4, 

+  4, 


+  '•'M.j-l  ^0  +  vM,j  10* 

/  J_  h__  ,  h_\ 

M-l.j+l  10  ‘  uM,j  40  vM,j  40 

,  1_  h_  h_> 

M-l  ,j-l  10  ‘  uM,j  40  ‘  vM,j  40' 

—  1  /I  h  \ 

V2,j  20  "  uM,j  40; 


■J  ;.SV  *  155 


' 
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I  /  It,  h  \ 

vM,j+2  ''  20  vM,j  W 


1  h 

+  *M,j-2  20  ~  vM,j 


h3  2 

120  ^  WM,j  WM,j+2  +  8wM,j+l  '  8wM , j - 1  +  wM,j-2^/ 

W 


r  i  1  ,  h  -| 

u  '  20  uM,j  40J 


At  the  point  i  =  2,  j  =  2,  ^  2  is 


h 


h 


*2,2  "  ^*3,2  ho  "  u2 ,2  10^  +  ^2 ,3  ho  "  v2,2  10^ 


+  *3,3  10  +  u2,2  40  +  v2 ,2  40^ 


+  *4,2  20  +  u2 ,2  40^ 


+  *2,4  20  +  v2 ,2  40^ 


h3  c2 

^720^7^  w2 ,2  ^w2,5  ‘  6w2,4  +  18w2,3  "  10w2,2^/ 


n  1_  h_  h_-| 

L  10  ‘  u2 ,2  40  "  2,2  40J 


At  the  point  i  =  M,  j  =  2,  ^  2  is 


(3-11) 


(3-12) 


*• 


4-U^  -  +  SH>'  t,MW  y '  W  ■ 


r  ;  s<r  ,S  t  .  J  ’  f  Jnr  9‘»  1 


m  S.SV  -F}  +  (UT  S,SU  'J{\  S.6'3  S.S* 


[m  S,SV  -F  S,SU  *  of  ■ 
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^M,2  =  ^M-1,2^10  +  UM,2  10^  +  ^,3  ^10  ‘  VM,2  lV 


+  ^-1,3  ("  10  "  ui,j  40  +  vi,j  40^ 


1  h 

+  ^M-2,2  20  ‘  ui,j  40^ 


+  ^M,4  20  +  vi ,j  40  ^ 


h3  wc2 


^120^r  WM,2  (wM,5  '  6wM,4  +  18wM,3  ~  10wM,2^ 


fi_L_+u  h _ v  J2 _ -| 

u  10  UM,2  40  VM ,2  40J 


At  the  points  i  =  3 ,4  . . .  M  -  1  ,  j  =  2 ,  if; .  ? 


becomes 


'h  ,2  =  ^1+1,2  ^0  '  ui  ,2  10^  +  ,3  *10  "  vi  ,2  10^ 


+  ’h-l ,2  <lo  +  ui ,2  10^ 


+  ^1+1 ,3  10  +  ui ,2  40  +  vi ,2  40  ^ 


,  /  J_  h__  h_ 

v1-l  ,3  v"  10  '  l ,2  40  vi  ,2  40 


1  h 

+  ’•'i-2,2  20  ‘  ui  ,2  40^ 


.  /  1_  h_ 

^i+2,2  K~  20  ui ,2  40 


(3-13) 


-  £,Mw8t  +  mw3  -  a,Mw)  S,MW  (7)(wr}  - 


>  z,^  +  (0T  S,ru  ■  0T:  S.f+f^  "  S.r 
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+  ,4  20  +  vi  ,2  40  ^ 


h3 

h  \  /  c 


-  (y2C))  (r  )  wi,2  (wi,5  ‘  6wi  ,4  +  18wi,3  ‘  10wi,2^/ 


n  1_  h _ -I 

L  '  20  "  i,2  40J 


(3-14) 


It  is  noted  that  in  the  derivation  of  equations  (3-10),  (3-11),  (3-12), 
(3-13)  and  (3-14)  for  equation  (2-28)  higher-order  approximation  is  used 
for  9w/9y  as  compared  with  the  inertia  terms. 

The  finite-difference  approximations  for  the  secondary  flow 
velocity  components  u,v  are: 

At  the  points  i  =  2,3  ...  M,  j  =  3,4  ...  N+  1 


=  (»j,j-2  -  +  8>N,j+1  ~  tyj  ,j+2^ 


ui,j  "  12h 


At  the  points  i  =  2,3  ...  M,  j  =  2 


(3-15) 


u 


i.2 


(*1,5  '  ,4  +  18l^i  ,3  ~  10^i  ,2^ 


W 


At  the  points  i  =  3 ,4  . . .  M  -  1  ,  j  =  2 ,3  . . .  N  +  1 


v .  . 
i  ,J 


=  ^i+2,j  '  8Vl,j  +  8Vl,j  ~  V2J) 

12h 


At  the  points  i  =  2,  j  =2,3  ...  N  +  1 


(Sf-I  rr-  )  ,(Ot-r)  ancrjti  is  to  nor jsv'neb  sri}  nr  aston  a  i 
at  *  rtMntxrmm*  nsbno- .arteM  (3S-S)  norJeups  not  (*r-£)  bne  r*<> 


.affrr  .J  s'rJinn'r  arlt  Urw  bansqmoo  6  tc«\w6  not 


vroft  vnsbnoosa  sriJ  not  ano'rJ6intxonqq£  aonane' .rb-at/n  t 


..  ..  -  S.S  “  >  into.,  «tJ  3A 


■  «7  '  « 
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2,j 


.  (»5.,1  ~  6^4  „j  +  1 8^3  „i  -  1(^2„i) 


-  1 2h 


(3-16) 


At  the  points  i  =  M,  j  =  2 ,3  . . .  N  +  1 


VM,j  = 


^M-3,j  '  ^M^.j  +  18Vl,j  ‘ 

12h 


The  finite-difference  expressions  for  the  boundary  conditions 
(equation  (2-31))  are: 

(1)  Along  the  lower  horizontal  surface, 
i  =  1  ,2  ...  M  +  1 ,  j  =  1 


wi  ,1  '  6i,l  '  ,j  =  ui  ,1  =  vi  ,1  "  0 


(2)  Along  the  horizontal  center  line, 
i  =  1 ,2  ...  M  +  1 ,  j  =  N 


wi,N  wi  ,N+2  ’  ei,N  ei,N+2  ’ 


^1  ,N  "  ,N+2  ’  .N-1  '  ’N.N+S 


’N.N+l  "  vi  ,N+1  0 


(3-17) 


(3)  Along  the  inner  wall, 


i  =  1 ,  j  =  1 ,2 


•  •  • 


N  +  1 


0  “  r.rv  “  r.tu  ‘  f.fr’ 


,9  t 


« rw  =  MW 


*  ’  =  |.M4'  '»t''  '  =  M 


' 


WV“  W‘  f\ 

. 


r  -u  .  .  s.r  =  t.f  -  r 
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w 


(4)  Along  the  outer  wall, 

1  =  M  +  1 ,  j  =  1 ,2  ...  N  +  1 


WM+1  ,j  =  ®M+1  ,j  =  Vi  ,j  =  UM+1  ,j  =  Vl  ,j  =  0 
3.2  Iterative  Method 

Point  successive-overrelaxation  method  [56]  is  used  to  solve 
a  set  of  finite-difference  equations  (3-6),  (3-9)  and  (3-7)  and  the 
boundary  conditions.  It  is  noted  that  the  two  momentum  equations  and 
the  energy  equation  are  not  coupled.  Consequently,  the  two  momentum 
equations  can  be  solved  independently.  The  iterative  procedure  employed 
is  similar  to  that  used  in  references  [51,53].  In  order  to  show  the 
procedure  clearly,  the  schematic  flow  chart  is  presented  below. 


Inner  Iteration 
for  Axial  Momentum 
Equation 


noH&islI  nsnnl 


■ 


,  fffiw  wel  0  ertf  gm> fA  (M 


bos  znorlsupe  roulnsmom  owJ  art*  JsrtJ  bslop  ?f  JI  .znorirbnob  'C**™0'1 

:  j 


®‘,9n9  9rtJ 
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where  eM,  and  eP  are  prescribed  error  bounds  for  w.  . ,  ip.  .  and  0.  . , 

not  i  jJ  I  $J  1  jJ 

respectively.  WK,  SK,  K,  and  TK  represent  the  number  of  sweeps  of  the 

inner  iterations  for  w.  ip-  that  for  outer  iteration  for  the  two 

i  *3  • 

momentum  equations  and  the  number  of  iterations  for  the  energy  equation, 
respecti vely .  TIMM,  TIMS,  TIMU,  and  TIME  represent  the  prescribed  numbers 
of  sweeps  for  the  corresponding  iterations  shown  in  the  chart. 


Inner  Iteration 
for  Secondary  Flow 
Stream  Function 


S'lSflW 


id  2  9  2  Temuri  dnses'rqn  X  -ne  «>i  «»:  «>'••' 
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The  error  is  inherent  in  the  numerical  solution.  The  errors 
due  to  finite-difference  approximations  for  the  differential  equations 
can  be  reduced  by  increasing  the  number  of  divisions  M  and  N  (see  Figure 
1).  A  typical  example  of  grid  size  effect  is  shown  in  Figure  2.  The 
number  of  divisions  M,  N  used  are  listed  in  Table  14.  The  round  off 
is  considered  to  be  insignificant  because  of  double  precision  (16  figures) 
and  a  point  iterative  method  used  in  this  study.  Numerical  experiment 
shows  that  the  round-off  error  is  quite  small  even  for  single  precision 
computation.  The  following  test  is  used  at  the  end  of  each  iteration  to 
determine  the  convergence  of  the  computed  function. 

Max|f.^  -  dn:1h/  Max|f.^|  <  e  (3-18) 

where  f.  .  may  represent  w.  . ,  ip.  . ,  e.  . ,  u.  .  and  v.  . ,  and  e  is  a 

1  )  J  1  )  J  1  )  J  1  5  J  1  j  J  1  )  J 

prescribed  error.  The  prescribed  errors  used  in  this  study  are  as  follows. 

_5 

e-.  =1  x  10  for  w.  . ,  iK  .  and  0  .  . 

I  1  $  J  •  I  j  J 

_5 

e9  =  5  x  10  for  u.  .  and  v.  . . 

The  above  error  bounds  give  sufficient  accuracy  up  to  a  reasonably  high 
Dean  number.  Specifically,  the  following  combinations  of  error  bounds  are 


e,  =  10'2,  10‘3,  10'4,  10"5,  10'6,  10'7,  and  10‘10 


examined. 


rugH  992)  H  bfi‘>  M  snofztvrb  >o  isdmun  srti  BnU*rt3iif  *d  fcsoubs  s 


Vto  bnuo'i  9rtT  .M  »rdsT  fit  b*f*H  *-b  beau  II  ,ft  zno^fvrb  to 


rtorri  rfdftnoe&a^  6  oi  qu  <?nu:>Dft  twbmw 
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e2  =  5  x  10'2,  5  x  10'3,  5  x  10"4,  5  x  10‘5,  5  x  10'6,  5  x  10"7  and  10'10. 

_5 

Comparison  of  the  flow  and  heat  transfer  results  using  e-j  =  1  x  10  , 

- 5  - 6  "6 

e2  =  5  x  10  and  e-j  =  1  x  10  ,  e2  =  5  x  10  ,  respectively,  shows  that 

the  difference  in  the  numerical  results  is  less  than  0.03%  even  in  high 

Dean  number  regime. 

For  the  number  of  divisions,  M,  N  used  in  this  study  (see 

Table  14),  the  flow  and  heat  transfer  results  obtained  by  two  different 

methods  (see  Section  4.1)  converge  to  nearly  the  same  numerical  value 

-5  -5 

with  e-j  =  1  x  10  and  =  5  x  10  as  can  be  seen  from  Tables  2  and 

2 

14.  For  example,  when  y  =  1,  M  =  32,  N  =  16  and  c  /r  =  0,  the  flow  and 
heat  transfer  results  obtained  by  two  different  methods  are  off  by  only 

-  0.0083%  and  0.0527%  for  flow  and  -  0.148%  and  0.105%  for  heat  transfer 

from  the  known  exact  values.  Furthermore,  when  y=l,M=32,N=16 

p  c 

and  c  /r  =  1.0  x  10  ,  the  flow  and  heat  transfer  results  using  two 
different  methods  are  off  by  0.220%  and  0.222%,  respectively.  These 
examples  confirm  the  accuracy  of  the  numerical  solutions. 

The  effect  of  grid  size  on  computing  time  is  shown  in  Figure  3. 

It  is  noted  that  this  effect  is  significant  and  the  number  of  divisions 

M,N  used  in  this  study  are  mainly  dictated  by  the  accuracy  of  the  numerical 

results.  A  relation  between  an  error  bound  and  the  computing  time  is  also 

of  considerable  interest  and  is  shown  in  Figure  4.  It  can  be  seen  that 

the  prescribed  error  used  in  this  study  is  quite  reasonable  from  the 
point  of  view  of  computing  time. 


992)  ybuds  nt  bszu  M  ,M  ,21-ofsrvfb  to  -ladmun  grid  -lo't 


' 


■'  ''4. 


jno  yd  Wo  916  !  o  i  J  )•  ■  ’  rb  owJ  yd  b9otbJdo  ti  <nar  istawl  J*9ri 


9i1  bnb  wort  9»  e  f  x  1 •  f  =•  3  •  .  0  >f- 


,,  )HT  .yfdv  toaqas'i  SSS.O  bne.XOSS.O  no  yte  dbaittan  Jnw' ‘  rb 


;UrU  nsea  sd  nfij  Ji  .*  S'luerl  nr  nworia  ar  bi  &  JaemJnr  sfds^sbra.noD  to 

,9fnrt  .tuqmoo  ‘to  W9‘rv  to  trro', 


40 


For  the  given  aspect  ratio  and  the  Prandtl  number,  the  numerical 
solution  starts  with  K  =  0  and  proceeds  gradually  toward  the  high  Dean 
number  regime.  The  convergence  of  the  iterative  method  is  closely  re¬ 
lated  to  the  suitable  choice  of  relaxation  factor  but  this  problem  will 
be  considered  in  the  next  section.  With  the  increase  of  the  Dean  number 
and  maintaining  the  same  prescribed  errors,  the  numerical  solution  is 
convergent  up  to  a  reasonably  high  Dean  number,  for  example,  for  a  curved 
square  channel;  however,  starting  at  a  certain  Dean  number  the  flow 
pattern  for  the  secondary  flow  changes  from  the  regular  two  vortices  to 
the  four  vortices  with  two  additional  small  vortices  located  near  the 
central  part  at  the  outer  wall.  The  results  of  numerical  experiment  at 
high  Dean  number  regime  are  shown  in  Figure  5  for  the  flow  characteristics. 
The  direction  of  the  numerical  solution  for  the  flow  is  indicated  by  the 
arrow  and  the  two  numbers  under  the  arrow  represent  the  number  of  inner 
iterations  for  the  axial  momentum  equation  and  the  secondary  flow  stream 
function,  respectively.  One  sees  that  a  different  combination  of  the 
number  of  inner  iterations  may  lead  to  a  different  flow  pattern.  In 
Figure  5,  two  different  flow  patterns  are  signified  by  a  circle  (2  vortices) 
and  a  triangle  (4  vortices),  respectively.  The  results  in  Figure  5 
seems  to  indicate  the  existence  of  a  pair  of  solutions  for  a  given  Dean 
number  in  high  Dean  number  regime.  Because  of  uncertainties,  the  flow 
and  heat  transfer  results  for  the  Dean  number  range  with  this  peculiar 
behavior  will  not  be  presented.  The  axial  velocity  and  temperature  pro¬ 
files  in  the  central  region  also  indicates  unusual  behavior  with  the  four 


.<  t  vT'  :-ofo  zt  boriin-  ovrJsisir  to  eone  i^vnoo  -  H\ 

,,  un  0  9  o'  lit  r:ni  srlJ  r  t'  orJ  a  1>«I  »«  <■  9  1  903  9(1 


"  ne9°  rieN 

.  ,  (zj*  1  *  •'•)  :/^  r^  t  o 

. 

. 

■ 


41 


vortices  for  the  secondary  flow. 

When  the  magnitude  of  the  secondary  velocity  components  exceeds 
2/h,  the  numerical  results  start  oscillation  for  the  solution  with  regu¬ 
lar  two  vortices  and  finally  the  solution  diverges.  The  above  situation 
occurs  at  a  certain  high  Dean  number.  One  notes  that  the  matrix  is 
diagonally  dominant  when  |u.  .|  or  |v.  .|  £  2/h.  On  the  other  hand,  the 
numerical  solution  with  four  vortices  is  still  convergent  even  when 
the  Dean  number  far  exceeds  the  value  indicated  above  for  the  flow  pattern 
with  two  vortices.  Eventually,  the  numerical  solution  with  four  vortices 
also  diverges  when  the  magnitude  of  the  secondary  velocity  components 
(u,v)  exceeds  2/h. 

Some  example  on  the  computing  time  required  may  be  of  interest. 

It  takes  about  138  minutes  by  IBM  360/67  to  obtain  a  complete  result  up  to 

c^/r  =  0.17  x  10^  for  flow  and  heat  transfer  with  y  =  1,  M  =  32,  N  =  16 
0 

and  Pr  =  0.71. 

3.3  Point  Successive-Overrelaxation  Method 

Since  point  iterative  methods  are  well  discussed  elsewhere 
[47,  52,  55,  56,  57],  the  details  of  point  successive-overrelaxation 
method  will  not  be  given  here.  To  improve  convergence  in  the  process 
of  iterations,  a  relaxation  factor  u>  is  used.  The  value  of  the  factor 
go  usually  lies  between  1  and  2.  A  question  naturally  arises  as  to  the 
optimum  value  of  the  factor  co  that  will  yield  a  maximum  rate  of  con¬ 
vergence.  Unfortunately,  no  general  method  is  available  for  the  evalu¬ 
ation  of  optimum  relaxation  factor  for  the  el  1 ipti c-type  partial  differ¬ 
ential  equations  having  the  non-linear  terms  encountered  in  the  present 


2bs90x9  a3n9fioqmoo  tfroorav  wrimoo"  art3  3o  «W  n9:  , 


,wci>  V16bn0092  9(13  103  2931310V 


,UB«  rti-rw  norJuToa  9rt3  io3  naWalftMO  iisU  tifua*  9ft:  ,dvS 


23n9noqmoD  tfroolav  vpfibnowz  9rt3  3o  bbuJrnesm  srt3  nsriw  «av»»>b  od 


. 


M  pu  3 f u29*t  9J9tqraob  6  nraPdo  o3  «\03£  MSI  yd  irtMhn  8£f  3uod*  mM  M 


-ortwsafe  b922U32fb  flaw  ais  zborUom  9vr36i93r  3nroq  aofliS 


229D0iq  9ri3  nf  939991 avnop  gvoiqraf  oT  .a<ari  nav'e  sd  3on  ffrw  bond am 
0ri3  03  2S  aaane  yl  h  n,  sn  nomnup  A  -S  una  f  ..sat  t*>  t»H  V  « 

.  *  I  , ,  r  nr\ 


problem.  However,  for  the  low  Dean  number  region,  one  would  expect  that 
the  method  described  in  [57]  may  be  applicable. 


In  this  study,  the  optimum  relaxation  factor  for  linear 
system  is  used  when  the  Dean  number  is  zero  or  sufficiently  small  as 
shown  in  Table  14.  By  using  w  =  1.75  ~  1.82,  considerable  computing 
time  is  saved.  For  the  Dean  number  ranging  from  small  to  intermediate 
values,  different  relaxation  factors  for  equations  (3-6),  (3-7)  and 
(3-9)  are  used  after  considerable  numerical  experiments.  For  example, 
for  a  curved  square  channel,  the  relaxation  factors  used  finally  for  the 
axial  momentum  equation,  secondary  flow  stream  function,  and  the  energy 
equation,  and  respectively  are  as  follows: 
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One  notes  that  at  high  Dean  number  region,  the  relaxation 
factor  oo=l  is  found  by  trial  and  error  to  be  the  best  value.  The  de¬ 
tails  of  solving  equation  (3-1)  and  equation  (3-8)  by  point  successive- 
overrelaxation  method  is  well  discussed  in  reference  [51]. 
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CHAPTER  IV 


FLOW  AND  HEAT  TRANSFER  RESULTS 

4.1  Flow  and  Heat  Transfer  Characteristics 

It  is  possible  to  obtain  the  expressions  for  the  product  of 
friction  factor  and  Reynolds  number,  fRe,  and  the  Nusselt  number,  Nu, 
by  considering  either  the  velocity  and  temperature  gradients,  respectively, 
along  the  channel  boundary  or  the  overall  force  and  energy  balances, 
respectively,  for  the  axial  length  dZ. 

(1)  The  Resistance  Coefficient 

Following  the  usual  definition  for  the  product  of  friction 
factor  and  Reynolds  number,  one  obtains 


tii  D 
(fRe)  =  — %  - 
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Considering  the  velocity  gradient  along  the  boundary,  one  obtains 
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Substituting  equation  (4-2)  into  equation  (4-1),  one  has 


(fRe)j 
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(4-3) 


A  force  balance  on  a  flow  element  with  axial  length  dZ  gives; 
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Consequently,  one  obtains 


(fRe) jj  =  ir 
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(2)  Nusselt  Number 

The  Nusselt  number  is  defined  by 


(4-5) 


(4-6) 


where  h  can  be  obtained  in  two  ways.  The  heat  flow  rate  through  the 


45 


channel  wall  with  axial  length  dZ  is 
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Thus,  the  Nusselt  number  can  be  expressed  as 


(Nu)j  = 


I (we) | 


(4-8) 


On  the  other  hand,  an  energy  balance  of  the  bulk  flow  with  axial  length 
dZ  yields , 
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Hence,  one  obtains 


(4-9) 


(4-10) 


The  single  or  double  integrations  required  for  the  evaluation 
of  the  mean  values  are  carried  out  by  using  Simpson's  rule.  For  the 
evaluation  of  the  first  derivatives  such  as  (|^)w  and  (-~)w  »  five-point 
formula  is  found  to  be  satisfactory.  The  foregoing  two  methods  of 
evaluating  (fRe)  and  Nu  affords  checking  of  the  accuracy  of  the  numerical 
resul ts . 

One  notes  that  the  parameter,  K,  defined  by  equation  (2-19)  can 
be  rewritten  as 


r2  1/2  _ 
K  =  (£-)  w 
c 


(4-11) 


and  will  be  used  in  the  following  sections. 

4.2  The  Effect  of  Dean  Number  on  Velocity  and  Temperature  Fields 

The  order  of  magnitude  analysis  shows  clearly  that  the  important 
characteristic  parameters  for  the  flow  and  heat  transfer  characteristics 
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are  Dean  number  and  Prandtl  number.  In  order  to  see  the  effect  of  Dean 
number  on  velocity  and  temperature  fields,  the  dimensionless  axial  velocity 
and  temperature  profiles  along  the  central  horizontal  axis  Y  =  0  and  the 
vertical  axis  X  =  0  are  plotted  for  a  curved  square  channel  y  =  1  with 
Pr  =  0.73  for  several  representative  Dean  numbers  in  Figures  6  and  7, 
respectively.  Qualitatively,  the  effect  of  centrifugal  force  on  the 
flow  and  temperature  fields  are  similar  to  the  effect  of  buoyancy  force 
[51]  for  a  given  geometrical  shape  of  channel  and  thermal  boundary  con¬ 
ditions  at  the  wall.  One  can  see  clearly  that  the  effect  of  the  centri¬ 
fugal  force  is  to  shift  the  location  of  the  maximum  value  toward  the  outer 
wall  and  decrease  the  maximum  value  itself  as  the  value  of  the  Dean  number 
K  increases.  It  is  also  seen  that  when  the  pressure  and  temperature 
gradients  in  the  axial  direction  are  kept  constant,  the  flow  rate  decreases 
and  both  the  velocity  and  temperature  profiles  become  flat  in  the  central 
region  of  the  channel  as  the  Dean  number  increases.  It  is  expected  that 
the  profiles  for  the  velocity  and  temperature  are  similar. 

Secondary  flow  streamlines  and  isothermals  for  a  curved  square 
channel  are  shown  in  Figure  8  for  K  =  51.9  and  Pr  =  0.73.  The  location 
of  the  center  of  circulation  is  of  interest  since  one  can  gain  the  general 
idea  about  the  secondary  flow  pattern  and  the  intensity  of  secondary  flow. 
For  a  curved  square  channel,  the  X-coordinate  of  the  center  of  circulation 
gradually  moves  from  X/a  =  0  toward  the  outer  wall  as  the  Dean  number 
increases  and  reaches  about  X/a  =  0.1  at  K  =  45.  With  further  increase 
of  the  Dean  number  the  center  of  circulation  tends  to  move  back  toward 
the  center  X/a  =  0.  It  is  found  that  with  Dean  number  at  about  125  the 
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center  of  circulation  returns  to  X/a  =0.  On  the  other  hand,  with  the 
increase  of  the  Dean  number  the  Y-coordinate  of  the  center  of  circu¬ 
lation  always  moves  toward  the  upper  or  lower  wall  indicating  the  in¬ 
crease  of  the  intensity  of  secondary  flow  near  the  upper  or  lower  wall. 

The  distributions  of  the  secondary  flow  velocity  components 
with  the  increase  of  the  Dean  number  are. also  of  considerable  interest 
and  these  are  shown  in  Figure  9.  One  also  sees  clearly  the  movement 
of  the  center  of  circulation  following  the  broken  lines  with  the  in¬ 
crease  of  the  Dean  number  in  this  figure.  One  notes  that  the  intensity 
of  the  secondary  flow  increases  as  the  value  of  the  Dean  number  increases. 

The  location  of  the  maximum  velocity  point  for  a  curved  square 
channel  is  also  of  interest.  As  the  Dean  number  increases,  the  location 
of  the  maximum  velocity  gradually  moves  toward  the  outer  wall.  At 
K  =  70,  the  maximum  velocity  is  located  at  around  X/a  =  .28  and  seems  to 
remain  there  with  further  increase  of  the  Dean  number.  One  may  add  that 
the  location  of  the  maximum  value  of  the  temperature  profile  shows  exactly 
the  same  trend  with  Pr  =  0.73. 

The  distributions  of  velocity  and  temperature  in  Figures  6  and 
7  suggest  that  boundary  layer  approximation  is  possible  for  both  velocity 
and  temperature  fields  when  the  Dean  number  is  greater  than  say  70  or 
preferably  125  on  the  conservative  side  (see  Figure  9).  This  observation 
is  important  since  it  confirms  the  applicability  of  the  boundary  layer 
approximation  for  the  high  Dean  number  region  [43].  For  the  low  Dean 
number  region,  the  boundary  layer  approximation  fails  and  the  present 
numerical  method  is  quite  effective. 
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One  can  gain  some  insight  into  the  flow  pattern  in  a  curved 
rectangular  channel  by  considering  the  distribution  of  centrifugal 
forces  and  pressure  in  a  cross-section.  One  notes  that  for  a  given 
radius  of  curvature,  the  centrifugal  force  is  proportional  to  the 
square  of  the  axial  velocity  at  a  given  point  and  acting  in  a  direction 
perpendicular  to  the  main  flow.  Consequently  the  fluid  in  the  central 
core  region  is  subjected  to  a  much  larger  centrifugal  force  than  the 
region  near  the  wall.  Due  to  the  centrifugal  force,  the  fluid  in  the 
central  core  will  be  pushed  toward  the  outer  wall  and  pressure  gradient 
arises  throughout  the  cross-section.  For  a  given  Y-coordinate,  the 
pressure  is  greatest  at  the  outer  wall  and  smallest  at  the  inner  wall. 

For  a  given  X-coordinate  to  the  right  of  the  center  of  circulation,  the 
pressure  is  greatest  at  Y  =  0  and  decreases  toward  the  upper  or  lower 
wall.  By  looking  at  the  secondary  flow  streamlines,  one  can  also  see 
the  distributions  of  pressure  gradients  throughout  the  cross-section. 

For  example,  the  strong  secondary  flow  toward  the  inner  wall  near  the 
upper  or  lower  wall  is  caused  by  large  pressure  drop.  Whereas  the  out¬ 
ward  flow  in  the  core  region  is  caused  by  centrifugal  forces.  It  is  well 
to  note  that  the  value  of  the  dimensionless  stream  function  at  the  center 
of  circulation  is  proportional  to  the  flow  rate  for  the  secondary  motion. 
When  the  circulation  in  a  cross-section  is  established,  the  main  flow 
pattern  will  be  distorted;  due  to  the  centrifugal  force,  the  location 
of  the  maximum  velocity  will  be  moved  toward  the  outer  wall.  For  a 
given  pressure  drop  in  the  axial  direction,  the  volume  flow  rate  de- 
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creases  because  of  the  secondary  motion.  Furthermore,  the  intensity  of 
the  secondary  motion  increases  with  the  increase  of  the  Dean  number.  The 
axial  velocity  distribution  with  the  increase  of  the  Dean  number  shown 
in  Figure  6  can  be  understood  readily  by  the  above  reasoning. 

4.3  The  Effect  of  Aspect  Ratio  on  Velocity  and  Temperature  Fields 

In  order  to  see  the  effect  of  aspect  ratio  on  flow  and  heat 
transfer  characteristics,  the  aspect  ratios  y  =  2,  5,  0.5  and  0.2  are 
considered  in  addition  to  a  curved  square  channel  y  =  1.  The  effect  of 
the  Dean  number  on  velocity  and  temperature  fields  in  a  curved  rectangular 
channel  with  various  aspect  ratios  is  generally  similar  to  that  for  a 
curved  square  channel . 

Consider  a  curved  rectangular  channel  with  long  side  horizontal 
first.  Dimensionless  axial  velocity  and  temperature  profiles  along  the 
central  axes  of  the  cross-section  of  the  channel  with  y  =  2  and  Pr  =  0.73 
are  shown  in  Figures  10  and  11,  respectively,  for  several  representative 
Dean  numbers.  Nearly  the  same  trend  exists  for  y  =  2  as  for  y  =  1. 

However,  for  a  given  value  of  Dean  number  the  locations  of  the  maximum 
velocity  and  the  maximum  value  of  the  temperature  profiles  move  closer 
to  the  outer  vertical  wall  as  compared  with  the  case  for  y  =  1 .  This  is 
due  to  the  distribution  of  centrifugal  forces.  One  notes  that  when  the 
Dean  number  exceeds  about  K  =  121.4,  the  temperature  profile  near  the 
central  region  of  the  channel  shows  slight  concave  upward  trend.  Secondary 
streamlines  and  isothermals  for  y  =  2,  Pr  =  0.73  and  K  =  58.8  are  shown  in 
Figure  12.  One  can  see  that  the  center  of  circulation  is  fairly  close  to 
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the  outer  wall  and  the  secondary  motion  is  very  weak  near  the  inner  wall. 
This  effect  is  also  reflected  in  the  isothermals.  The  distributions  of 
the  secondary  velocity  components  (u,v)  are  shown  in  Figure  13  with 
y  =  2  for  several  Dean  numbers.  The  loci  of  the  center  of  circulation 
with  the  increase  of  the  Dean  number  are  also  shown  in  the  same  figure. 
The  increase  of  the  intensity  of  secondary  flow  with  the  increase  of  the 
Dean  number  is  quite  evident. 

For  the  case  y  =  5  and  Pr  =  0.73,  graphical  results  for  axial 
velocity  and  temperature  profiles  with  K  as  a  parameter,  and  secondary 
flow  streamlines  and  isothermals  for  K  =  88.1  are  shown  in  Figures  14 
to  16.  Once  again  one  sees  that  secondary  motion  is  quite  strong  near 
the  outer  wall  and  rather  weak  near  the  inner  wall.  This  effect  is  again 
reflected  in  the  isothermals.  The  locations  of  the  maximum  velocity  and 
the  maximum  value  of  the  temperature  profile  are  much  closer  to  the  outer 
wall  for  y  =  5  than  y  =  2  for  the  same  Dean  number. 

The  effect  of  the  aspect  ratio  when  the  long  side  is  vertical 
will  be  examined  next.  The  velocity  and  temperature  distributions  with 
K  as  a  parameter  along  the  central  axes  of  a  curved  rectangular  channel 
with  y  =  0-5  and  Pr  =  0.73  are  shown  in  Figures  17  and  18,  respectively. 
Secondary  flow  streamlines  and  isothermals  for  the  aspect  ratio  y  =  0.5, 
Pr  =  0.73  and  K  =  103,4  are  presented  in  Figure  19.  The  secondary  flow 
patterns  for  several  representative  Dean  numbers  and  the  loci  of  the 
center  of  circulation  are  shown  in  Figure  20.  One  notes  particularly 
that  with  the  increase  of  the  Dean  number,  the  center  of  the  circulation 
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tends  to  move  toward  the  upper  or  lower  horizontal  wall;  however,  the 
center  tends  to  move  slightly  toward  the  outer  wall  with  the  increase  of 
the  Dean  number  and  after  reaching  about  K  =  70,  tends  to  move  toward 
the  inner  wall  with  further  increase  of  the  Dean  number.  It  is  expected 
that  for  the  high  Dean  number  regime,  the  intensity  of  the  secondary  flow 
will  be  quite  strong  near  the  upper  or  lower  wall  as  compared  with  the 
central  core  region. 

Similar  graphical  results  for  the  case  y  =  0.2  are  shown  in 
Figures  21  to  23.  For  this  case  with  the  increase  of  the  Dean  number 
the  center  of  circulation  moves  toward  the  upper  or  lower  wall  and  in¬ 
tensive  secondary  flow  exists  only  near  the  upper  or  lower  wall.  Con¬ 
sequently,  for  the  low  Dean  number  regime,  secondary  motion  has  less 
effect  on  the  velocity  and  temperature  distributions.  Figure  23  shows 
that  at  K  =  79.3,  the  center  of  the  circulation  is  located  quite  close 
to  the  central  vertical  axis  of  the  cross-section  but  slightly  on  the 
outer  wall  side.  Although  the  Dean  number  is  fairly  high,  the  secondary 
flow  is  weak  throughout  the  central  region.  It  is  not  difficult  to  see 
the  general  trend  for  the  velocity  and  temperature  fields  with  further 
decrease  of  the  aspect  ratio,  namely  as  y  +  0.  However,  one  notes  that 
as  y  +  0,  the  problem  leads  to  the  instability  problem  discussed  in 
reference  [45]. 
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4.4  The  Effect  of  Prandtl  Number  on  Temperature  Field 

Observation  of  the  equations  of  motion  (2-28)  and  (2-29)  shows 
that  Prandtl  number  has  no  effect  on  the  flow  field.  In  order  to  under¬ 
stand  the  effect  of  Prandtl  number  on  temperature  field  equations  (2-29) 
and  (2-30)  will  be  rewritten  as  follows. 

Axial  Momentum  Equation: 


u 


3w  ,  3w 
3x  v  3y 


V2w 


4 


Energy  Equation: 


=  v2e-w 

First,  one  notes  immediately  the  similarity  between  the  above  two  equations. 
In  fact  if  Pr  =  1 ,  the  axial  velocity  distribution  is  similar  to  the 
temperature  distribution  for  a  given  Dean  number.  Furthermore,  the  effect 
of  inertia  terms  in  the  axial  momentum  equation  increases  with  the  increase 
of  the  Dean  number.  Consequently,  the  effect  of  Prandtl  number  on  the 
left-hand  side  of  the  energy  equation  is  similar  to  the  effect  of  Dean 
number  on  the  left-hand  side  of  the  momentum  equation.  Noting  the  effect 
of  the  Dean  number  on  the  axial  velocity  distributions,  one  can  immediately 
see  the  effect  of  Prandtl  number  on  the  temperature  distributions.  The 
effect  of  Prandtl  number  can  be  seen  from  Figures  24  to  26  where  temperature 
distributions  through  the  central  axes  of  a  curved  square  channel  y  =  1 
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2  3 

with  K  as  a  parameter  are  shown  for  Pr=  0.1,  10  and  10  ,  respectively. 
The  Prandtl  number  effect  explained  above  is  clearly  seen. 

4.5  Flow  Resistance 

The  overall  flow  characteristics  will  be  considered  next.  The 
ratio  of  the  average  velocity  of  a  curved  rectangular  channel  to  that  of 
a  straight  channel  is  plotted  against  Dean  number  for  the  aspect  ratios 
y  =  0.2,  0.5,  1,  2  and  5  in  Figure  27.  One  sees  that  for  a  given  axial 
pressure  gradient,  the  effect  of  Dean  number  on  average  axial  velocity 
is  greatest  for  the  curved  square  channel  up  to  K  ~  100.  One  also  notes 
that  as  the  aspect  ratio  y  approaches  zero  or  infinity,  the  Dean  number 
effect  decreases.  One  can  see  from  the  figure  that  the  curves  for  y  =  1 
and  2  cross  at  about  K  =  100  and  similarly  the  curves  for  y  =  5  and  0.5 
also  cross  at  some  value  of  K. 

The  ratio  of  the  products  of  friction  factor  and  Reynolds  number 
(fRe)  between  a  curved  rectangular  channel  and  a  straight  rectangular 
channel  is  plotted  against  Dean  number  in  Figure  28  for  the  aspect  ratios 
y  =  0.2,  0.5,  1,  2  and  5.  Of  course  the  product  (fRe)  and  the  average 
velocity  are  related  to  each  other  as  seen  from  equation  (4-5).  It  is 
of  interest  to  note  that  the  aspect  ratios  y  =  0.5  and  2  represent  the 
same  cross-sectional  area.  Similarly,  the  cross-sectional  area  is 
identical  for  the  aspect  ratios  y  =  0.2  and  5.  Comparison  of  the  curves 
for  y  =  2  and  0.5  in  Figure  28  shows  that  the  Dean  number  effect  is  much 
stronger  for  a  channel  with  larger  aspect  ratio  after  reaching  a  certain 
value  of  K.  Similar  remark  applies  for  the  cases  y  =  5  and  0.2  with  the 
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same  cross-sectional  area.  This  fact  may  be  of  interest  in  design.  When 
K  is  small,  the  centrifugal  force  effect  is  small  and  the  inertia  terms 
may  be  negligible  as  compared  with  the  viscous  terms;  consequently,  the 
way  of  placing  long  side  horizontal  or  vertical  has  negligible  effect  on 
flow  resistance.  In  contrast,  for  the  high  Dean  number  region,  the 
effects  of  the  centrifugal  force  term  and  the  inertia  terms  are  signi¬ 
ficant.  The  different  effect  of  the  aspect  ratio  with  the  same  cross- 
sectional  area  may  be  explained  from  the  distribution  of  the  centrifugal 
forces  for  the  high  Dean  number  regime. 

The  results  for  the  various  aspect  ratios  presented  in  Figure 
28  show  that  for  the  high  Dean  number  regime,  the  value  of  (fRe)/(fRe)Q 
changes  Ti nearly  with  Dean  number  K;  consequently,  one  may  write 
(fRe)/(fRe)QCO Km  where  m  is  a  function  of  the  aspect  ratio.  For  example, 

(fRe)/(fRe)0  *  0.225  K0-39  for  y  =  1  and  102  <  K  <  1.5  x  103 

Comparison  of  the  result  from  the  present  analysis  with  the  results 
available  in  the  literature  is  of  considerable  interest  and  is  shown  in 
Figure  29.  One  can  see  that  the  present  study  covers  the  Dean  number 
ranging  from  small  to  a  reasonably  high  region  where  no  other  work  is 
available  in  the  literature.  One  numerical  data  at  a  quite  high  Dean 
number  is  also  plotted  in  Figure  29  for  comparison  and  shows  a  very  good 
agreement  with  Ludwieg's  experimental  data  [38].  Mori  and  Uchida's 
analysis  [43]  using  boundary  layer  approximation  for  the  high  Dean  number 
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regime  agrees  well  with  Ludwieg's  experimental  data  up  to  a  certain  Dean 
number.  Beyond  that  the  experimental  data  by  Ludwieg  are  suspected  to 
be  in  the  turbulent  region.  The  two  curves  given  by  Mori  and  Uchida 
represent  the  first  and  second  approximations. 

Observation  of  the  result  from  the  present  numerical  analysis 
and  the  results  from  Ludwieg  [38]  and  Mori  and  Uchida  [43]  shows  clearly 
that  a  reasonable  estimate  can  be  made  for  the  flow  resistance  with  the 
Dean  number  ranging  from  150  to  1000  where  currently  accurate  analytical 
solution  is  not  available. 

4.6  Heat  Transfer 

The  overall  heat  transfer  characteristics  will  be  examined  next. 
In  order  to  see  the  effect  of  Dean  number,  the  ratio  of  the  average 
temperatures  F/(0")q  between  a  curved  channel  and  a  straight  channel  is 
plotted  against  K  in  Figure  30  for  y  =  0.2,  0.5,  1,  2  and  5  with  Pr  = 

0.73.  One  notes  that  the  Dean  number  effect  is  quite  similar  for  0/ ( 0~) q 
as  for  w/(w)Q 

In  order  to  see  the  Prandtl  number  effect  on  ?/(¥)q,  0/(?)q  is 
plotted  against  Dean  number  K  for  a  curved  square  channel  y  =  1  in 
Figure  31  for  Pr  =  Q.l,  0.73,  1,  10,  10^  and  10^.  One  can  see  that 
Prandtl  number  effect  is  considerable. 

The  graphical  results  for  the  ratio  of  the  Nusselt  number  of  a 
curved  channel  to  that  of  a  straight  channel  versus  Dean  number  are 
shown  in  Figure  32  for  the  aspect  ratios  y  =  0.2,  0.5,  1,  2  and  5  with 
Pr  =  0.73.  As  can  be  expected,  heat  transfer  result  is  similar  to  flow 
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resistance  result  discussed  earlier.  As  noted  earlier,  the  effect  of 
Prandtl  number  on  temperature  field  is  considerable.  This  is  also  re¬ 
flected  in  heat  transfer  result  shown  in  Figure  33  where  Nu/(Nu)q  is 
plotted  against  K  for  a  curved  square  channel  with  Pr  =  0,  0.1,  0.71, 

0.73,  1,  10,  102,  102  and  10^.  The  results  in  Figure  33  show  clearly 
that  an  asymptotic  line  exists  for  the  heat  transfer  result,  too,  in 
the  high  Dean  number  region  with  a  given  Prandtl  number;  consequently, 
in  the  high  Dean  number  regime,  the  following  expression  may  be  used. 

Nu/(Nu)q  *  0.182  K1/2  Pr1/4  for  1.0  <  Pr  <  104 

and  Nu/(Nu)0  >_  1 .5 

It  is  of  interest  to  note  that,  for  example,  at  Nu/(Nu)q  =  1.6  the  distance 
between  two  neighboring  curves  decreases  slightly  as  Pr  increases  by  the 
same  factor. 

Again  it  is  of  considerable  interest  to  compare  the  heat  transfer 
result  from  the  present  analysis  with  the  result  from  Mori  and  Uchida 
[43]  for  a  curved  square  channel  y  =  1.  Mori  and  Uchida  [43]  show  the 
first  and  second  approximations  for  the  ratio  of  Nusselt  numbers  for 
Pr  =  0.71  and  °°.  The  comparison  of  the  heat  transfer  results  is  shown 
in  Figure  34.  A  data  from  numerical  solution  at  K  =  460  is  also  shown 
in  the  figure  for  comparison.  Concentrating  on  the  case  Pr  =  0.71,  one 
sees  that  the  present  result  is  very  reasonable  up  to  a  fairly  high  Dean 
number  and  suggests  clearly  that  a  reasonable  estimate  can  be  made  for 
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the  Dean  number  ranging  from  K  =  150  to  1000.  The  estimated  line  is  shown 
in  Figure  34  as  a  broken  line.  One  notes  the  asymptotic  nature  of  the 
heat  transfer  result  for  the  high  Dean  number  regime. 

Mori  and  Uchida's  work  [43]  shows  that  an  asymptotic  value  exists 
for  the  ratio  Nu/(Nu)q  as  Pr  however,  this  behavior  of  the  effect 
of  Prandtl  number  is  difficult  to  comprehend  in  view  of  the  energy  equation 
(2-30).  The  numerical  result  for  heat  transfer  from  this  study  shows 
that  an  asymptotic  line  for  Pr  =  “will  not  be  reached  at  least  up  to  the 
range  (Pr  =  10^)  studied  in  this  work.  Furthermore,  in  reference  [30] 
it  is  stated  that  the  Nusselt  number  ratio  approaches  the  asymptotic 
value  with  the  increasing  Prandtl  number  for  a  similar  problem  in  curved 
pipes.  Based  on  the  result  of  this  work,  it  appears  that  an  asymptotic 

4 

value  does  not  exist  at  least  up  to  Pr  =  10  .  On  the  other  hand,  an 
asymptotic  value  does  exist  as  Pr  ->  0. 

It  is  known  that  as  the  Dean  number  increases,  the  effect  of 
the  convective  terms  in  the  energy  equation  (2-30)  increases  for  a  given 
Prandtl  number.  As  the  Prandtl  number  increases,  the  effect  of  the  con¬ 
vective  terms  also  increases  for  a  given  Dean  number.  One  can  see  that 
the  Nusselt  number  increases  with  the  increase  of  the  Prandtl  number 
even  if  the  Dean  number  is  held  constant.  The  effect  of  the  Prandtl 
number  is  equivalent  to  the  effect  of  the  Dean  number.  For  example,  at 
Nu/(Nu)q  =  1.6,  as  the  Prandtl  number  increases  from  10  to  10  ,  the  Dean 
number  changes  from  K  =  20  to  K  =  2. 

In  addition  to  the  graphical  results  shown  above,  the  numerical 
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results  are  tabulated  in  Tables  1  to  15.  One  may  be  particularly  inter 
ested  in  the  numerical  results  based  on  two  different  methods  of  ob¬ 
taining  flow  and  heat  transfer  results  discussed  in  Section  4.1.  It 
can  be  seen  from  Tables  2  and  14  that  the  results  from  the  two  methods 
are  off  by  only  0.4%  at  most  for  all  the  cases  considered.  This  con¬ 
firms  the  accuracy  of  the  numerical  solutions.  The  numerical  values 
given  in  Tables  3  to  13  are  obtained  by  using  equations  (4-3)  and  (4-8) 
Comparisons  of  equation  (4-3)  with  equation  (4-5)  and  equation  (4-8) 
with  equation  (4-10)  show  that 


(§*■)  =  1  and 

3n  W 


36 1  w 

3n  w  "  4 


The  numerical  results  in  Table  15  also  confirms  the  accuracy  of  the 


numerical  solutions. 


(hjrwttw,  sd  *6m  snO  .at  O.t  T  «f*T  nf  fc^efud^  «?»  •trim*' 


. 


CHAPTER  V 


CONCLUDING  REMARKS 

1.  Numerical  solution  by  point  successive-overrelaxation  is  ob¬ 
tained  for  the  fully  developed  laminar  forced  convection  problem  in 
curved  rectangular  channels  with  various  aspect  ratios  under  the  thermal 
boundary  conditions  of  constant  heat  flux  per  unit  axial  length  with  a 
uniform  wall  temperature  for  the  conducting  surfaces  at  each  cross- 
section.  The  limitation  of  the  numerical  method  is  encountered  at  a 
reasonably  high  Dean  number.  The  difficulty  seems  to  come  from  the  non¬ 
linear  terms  in  the  governing  equations.  In  spite  of  this  difficulty, 
the  numerical  method  has  definite  advantage  over  the  perturbation  method. 

2.  Forced  convection  problems  in  rotating,  bent,  and  heated  hori¬ 
zontal  tubes  are  known  to  be  similar  and  characterized  by  the  secondary 
flow  acting  in  a  cross-section  normal  to  the  main  flow.  Buoyancy  force 
effects  on  forced-convection  flow  and  heat  transfer  in  horizontal 
rectangular  channels  with  various  aspect  ratios  are  studied  in  references 
[51,52].  The  aspect  ratios  considered  in  this  study  are  identical  to 
those  in  references  [51,52]  and  the  direct  comparison  is  possible  since 
the  thermal  boundary  conditions  are  also  identical.  Comparison  of  the 
results  shows  that  a  striking  similarity  exists  between  the  centrifugal 
force  effect  and  the  buoyancy  force  effect  on  the  flow  and  heat  transfer 
characteristics  within  the  range  studied  by  the  references  [51,52], 
namely,  Pr  =  0.71  and  7,  (fRe)/(fRe)Q  <  1.18  and  Nu/(Nu)Q  <  1.3.  One  is 
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impressed  with  this  similarity  in  view  of  the  fact  that  for  the  forced 
convection  problem  with  buoyancy  effect  the  momentum  equations  for  the 
main  flow  and  the  secondary  flow  are  coupled  with  the  energy  equation 
whereas  for  the  forced  convection  problem  in  a  curved  channel,  the  two 
momentum  equations  and  the  energy  equation  are  uncoupled.  It  is  noted 
that  the  effect  of  the  secondary  flow  caused  by  Coriolis  force  on  forced 
convective  heat  transfer  in  a  rotating  channel  is  recently  studied  in 
reference  [58]. 

3.  By  using  co  =  1  for  the  relaxation  factor  at  the  high  Dean 
number  region  and  at  the  sacrifice  of  the  computing  time,  the  range 
of  the  numerical  solution  in  terms  of  the  ratios  (fRe)/(fRe)Q  and 
Nu/(Nu)q  is  further  extended  beyond  the  range  studied  by  the  reference 
[51].  One  may  note  that  the  range  of  the  numerical  solution  in  reference 
[51]  may  be  extended  by  using  the  relaxation  factor  similar  to  the  one 
employed  in  this  study.  The  relaxation  factor  used  in  reference  [51]  is 
an  optimum  value  based  on  linear  system.  For  the  low  Dean  number  region, 
the  relaxation  factor,  1.5  ~  <  co  <  1.82,  based  on  linear  equation  may 

be  used.  For  a  curved  rectangular  channel  with  y  =  2  the  numerical  solu¬ 
tion  is  obtained  up  to  the  ratios  (fRe)/(fRe)^  =  1.8  and  Nu/(Nu)q  =  2.2, 
respectively.  One  sees  clearly  that  the  range  of  applicability  of  the 
present  numerical  solution  far  exceeds  that  of  the  perturbation  method 
[34]. 

4.  The  result  of  the  present  numerical  analysis  shows  clearly 
that  boundary-1 ayer  approximation  around  the  wall  for  the  high  Dean 
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number  regime  is  possible  and  it  clarifies  the  assumption  that  at  high 
Dean  number  a  flow  field  may  be  divided  into  a  core  region  occupying 
most  of  the  cross-section  where  the  viscous  stress  may  be  neglected  and 
the  boundary  layer  region  adjacent  to  the  wall  where  the  viscous  stress 
cannot  be  neglected.  The  numerical  method  is  effective  for  the  flow 
regime  with  Dean  number  ranging  from  small  to  a  reasonably  high  value 
and  complements  the  boundary  layer  technique  [43]. 

5.  Because  of  geometrical  symmetry,  two  vortices  with  opposite 
sense  are  formed  for  the  secondary  flow.  For  a  given  aspect  ratio  of 
the  curved  rectangular  channel,  the  intensity  of  secondary  motion  in¬ 
creases  with  the  increase  of  the  Dean  number.  The  secondary  motion  is 
independent  of  the  Prandtl  number.  For  a  given  Prandtl  number,  the 
effects  of  Dean  number  on  axial  velocity  and  temperature  profiles  are 
similar  as  can  be  seen  from  the  axial  momentum  equation  and  the  energy 
equation.  As  the  Dean  number  increases,  the  center  of  circulation 
always  tends  to  move  toward  the  upper  or  lower  wall.  The  center  of 
circulation  gradually  moves  toward  the  outer  wall  with  the  increase  of 
the  Dean  number  and  after  reaching  a  certain  Dean  number  tends  to  move 
toward  the  inner  wall.  The  location  of  the  center  of  circulation  is  of 
special  interest  since  this  gives  some  insight  into  the  secondary  flow 
pattern.  For  a  given  Dean  number,  the  effect  of  aspect  ratio  is  con¬ 
siderable.  Even  for  the  same  cross-sectional  area,  the  secondary  flow 
pattern  changes  considerably  depending  on  whether  the  long  side  is  hori¬ 
zontal  or  vertical.  This  fact  is  of  considerable  interest  in  design. 
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The  numerical  results  show  an  asymptotic  nature  for  both  (fRe)/(fRe)Q 
and  Nu/(Nu)q  in  the  high  Dean  number  regime. 

6.  The  effect  of  Prandtl  number  on  heat  transfer  result  is  signi¬ 
ficant.  The  role  of  the  Prandtl  number  can  be  seen  clearly  from  the 
energy  equation.  Considering  the  convective  terms  in  the  energy  equation, 
one  sees  that  the  effect  of  the  Prandtl  number  is  equivalent  to  the 
effect  of  Dean  number.  It  is  pointed  out  in  references  [30,43]  that  the 
Nusselt  number  ratio  Nu/(Nu)q  approaches  the  asymptotic  value  as  the 
Prandtl  number  approaches  infinity.  In  view  of  the  above  reasoning  and 
the  numerical  results  obtained  in  this  study  at  high  Prandtl  number, 

it  appears  that  the  heat  transfer  results  for  Pr  =  °°  shown  in  references 
[30,43]  are  incorrect.  In  reference  [30],  the  heat  transfer  results  ob¬ 
tained  from  boundary-layer  approximation  for  Pr  =  »  are  shown  to  agree 
with  the  experimental  results  obtained  by  Seban  and  McLaughlin  [28]  for 
Pr  ~  400.  It  is  believed  that  the  use  of  the  analytical  results  for 
Pr  =  oo  presented  in  references  [30,43]  will  lead  to  serious  error  in 
heat  transfer  prediction  when  Pr  >  400.  On  the  other  hand,  as  Pr  ap¬ 
proaches  zero,  the  heat  transfer  result  approaches  the  asymptotic  value 
as  shown  in  Figures  33  and  34.  One  may  note  that  the  boundary-layer 
approximation  cannot  be  used  as  Pr  -►  0. 

7.  Comparison  of  the  results  from  this  analysis  and  the  results 
available  in  the  literature  [38,43]  for  a  curved  square  channel  shows 
that  a  reasonable  estimate  can  be  made  for  the  flow  and  heat  transfer 
results  for  the  Dean  number  ranging  from  150  to  1000  where  currently 
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accurate  solutions  are  not  available. 

8.  The  order  of  magnitude  analysis  establishes  the  conditions 
under  which  the  present  analysis  may  apply.  Specifically,  it  shows 
under  what  conditions  the  buoyancy  force  effect  and  the  Coriolis  force 
effect  may  be  neglected  for  the  forced  convection  problem  in  a  curved 
rectangular  channel . 

9.  The  result  of  the  present  analysis  confirms  the  known  fact 
that  the  local  heat  transfer  coefficient  is  higher  at  the  outer  wall 
of  the  curved  channel  than  at  the  inner  wall. 
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APPENDIX 


FORTRAN  PROGRAM 


A  «*»  IN  COL.  6  INDICATES  THAT  THIS 
LINE  IS  A  PART  OF  THE  PREVIOUS  CARD 


DOUBLE  PRECISION  ERR , FM , FN , H , H2 , H4 , H3 » H4 , H6 , H8 , H9 

*  , H24 , H  72  » 

1W(81,43),T(81,43),S(81,43),U(81,43),V(81,43), 

2G AMMA , P R , CON , D I , RE  F , WM , TM , TW , FRE , T  3 , NU , ER 1 , VTE  S T , CGW2 J 

*  ,CGWMJ, 

3WA  f WNfWD  »WER  »TA, TN, TD,HA,TER  »RCTt S A , S N , SER , VN ,  UN , UD ,  Vu 

*  » 

4SH<81  ,43)  t TGI  81 ,43)  , SUM, SUT , SUM  I , SUT 1 , S H M , T GM , HB , VE R , S 
^  D , 

5Bl»B2»B3»B4,85,B6»B7f88,89,B10,Bll,812»B13»DK»RC*NUE,F 

*  P.EM, 

6  W  1 ,  W  2  ,  W  3  ,  U 1  ,  U2  ,  U3  *  U4  ,  U5  ,  V 1  ,  V2  ,  V3  ,  V4  ,  V 5  ,  OME  M, OMES  »  OMEfc 

*  » 

7H10 »H23»H40,H80,H120»H144»H240 ,HI 44  C,B40»B80 
INT  EGER  TIMM,T IMS, TIMU, TIME, WK,TK,SK 
COMMON  M 1 , N 3 

C  M  IS  THE  NUMBER  OF  DIVISIONS  IN  THE  X  DIRECTION  BLTWEE 

*  N  0  AND  A. 

C  N  IS  THE  NUMBER  OF  DIVISIONS  IN  THE  Y  DIRECTION  BETWEE 

*  NO  AND  312. 

C  GAMMA  IS  THE  ASPECT  RATIO  OF  A  RECTANGULAR  CHANNEL 

C  OME M, OMES  AND  OMEE  ARE  THE  RELAXATION  FACTORS 

C  FOR  Z_MOMENTUM, STREAM  FUNCTION  AND  THF  ENERGYEQU ATI  ON, 

C  RESPFECT IVLY. 

C  ERR  IS  THF  PRESCRIBED  ALLOWABLE  ERROR 

C  TIMM,  TIMS,  TIMU  AND  TIME  ARE  THE  PRESCRIBED 

C  NUMBERS  OF  SWEEPS  FOR  Z_MOMENTUM,  STREAM  FUNCTION, 

C  THE  VELOCITY  COMPONENTS  OF  THE  SECONDARY  FLOW  AND  THE 

C  ENERGY  EQUATION,  RESPECTIVELY. 

READ! 5,45 )  M, N, OMFM, OMES , OMEE , ERR  ,T I  MM, T I  MS , T I  MU, T I  ME 

M I=M+ I 
M I  =  M-  I 
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M j=m-2 
MK=  M-3 
N 1  =  N+  1 
N2=N+2 
N  3=N  +  3 
N I  =  N- 1 

FM=D3LE( FLOAT { M) ) 

FN=DBLE ( FLOAT ( N  >  ) 

C  U*V  AND  W  ARE  THE  DIMENSIONLESS  VELOCITY  COMPORNENTS  I 

*  N  THE  X,  Y  A 
Z  DIRECTIONS. 

T  IS  THE  DIMENSIONLESS  TEMPERATURE 
S  IS  THE  STREAM  FUNCTION 

I  AND  J  ARE  THE  SPACE  SUBSCRIPTS  OF  GRID  POINTS  IN  THE 
X  AND  Y  DIRECTIONS. 

DO  1  I  =  1 ,  M 1 
DO  1  J  =  1 ,  N  3 
W( I , J)=ODO 
T ( I , J)=ODO 
S< I t J)=0Q0 
U( If J)=ODO 
1  V( I , J)=0D0 

CALL  DATA  (W,T,S,U,V) 
bS5  CONTINUE 

K  IS  THE  NUMBER  OF  SWEEPS  OF  U  AND  V. 

PR  IS  THE  PRANDTL  NUMBER. 

CON  IS  THE  RATIO  OF  THE  SQUARE  OF  PRESSURE  GADIENT,C, 

*  TO 

THE  RADI OUS  OF  THE  CURV ETUR E , R C . 

OF  CURVETUREtRC. 

H  IS  THE  DIMENSIONLESS  GRID  SPACING  IN  THE  X  AND  Y  DIP 

*  EC T IONS. 

K=  1 

R  EA  DI 5 , 46 )  GAMMA, PR, CON 
I F ( GAMMA )  600,1000,600 
600  B 1=4D- 1 
62=10-1 
B3=-B2 
B5=-5D-2 
B40  =  2  50-3 
880=5D-1*B40 
H= ( GAMMA+1D0 ) /( 2D0*FM ) 

H2=H*H 
H3=H**3D0 
H4  =  H23!t25D-2 
H8=H/8D0 
H10=H*82 
H23=H3*2D-1 /3D0 
H40=H^340 
H120=H3/12D1 

WRITE! 6, 49)  M , N,  E RR ♦ G AMMA , PR , CON 
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ITERATION  PROCEDURE 


THE  START  OF  ITERATIONS 

m  IS  THE  NUMBER  OF  SWEEPS  OF  MOMENTUM  EQUATION  FOR  CO 

*  M  PUT  AT  ION  OF 

C  VELOCITY  PROFILE  WITH  RESPECT  TO  THE  AX  I AL , Z , D I RECT I  ON 

*  . 

24  WK=  1 

WR  ITE ( 6t  64 ) 

7  D I =  ODO 
REF -CDO 
DO  4  1  =  2,  M 
I  1=  I  + 1 
11=1-1 
DO  4  J  =  2  » N 1 
J  1=J  +  1 
JI=J-1 

U1=U< I ,  J )*H8 
V 1  =  V (  I  ,  J  ) *H8 

C  THE  VELOCITY  WITH  RESPECT  TO  THE  AXIAL  DIRECTION,  WA, 

*  I S 

WA  =  W {  II  ,  J)*<  250-2-U1 )+W<II  , J)*( 25D-2+U1 ) 

1  +W(  I  ,  J1  )*i 25D-2-V1 )+W(  I  ,  J  I  )*<  2  5D-2+V1 ) 

1  +H2 

WN  =  0MEM5MWA-w(  I,J)  )+W(If  J) 

I F ( J .LT.N.OR. J.GT.N)  GO  TO  5 
W( I  *  N  2 ) =WN 

5  WD=DABS ( WN-wI  I , J  )  ) 

C  D1  IS  THE  MAXIMUM  DIFFERENCE  BETWEEN  THE  VALUE  JUST  OB 

*  TAINED  AND  P 

C  VALUE 

DI=DMAX1( DI ,WD  ) 

C  REF  IS  THE  MAXIMUM  ABSOLUTE  VALUE 

REF=DMAX1 { REF,WN ) 

4  W ( I , J)=WN 
WER=DI /REF 

WRITE(6,51)  WK , WE  R , OMEM 
IFIWER.LE. ERR. GR.WK.GE. TIMM)  GO  TO  6 
WK=WK+ 1 
GO  TO  7 

6  CONTINUE 

C  SK  IS  THE  NUMBER  OF  SWEEPS  OF  THE  VORTICITY  (  STREAM  F 

*  UNCTION  ) 

C  EQUATION  USING  THE  OBTAINED  W  TO  COMPUTE  STREAM  FUNCTI 

*  ON 
$K=  1 

WRITE (6, 65) 

16  D I =ODO 
REF=ODO 
DO  12  1=3, MI 
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11=1+1 
12=1+2 
11=1-1 
I  J= 1-2 
DO  12  J  =  3  T  N 
J  1  =  J  +  1 
J  2  =  J  +  2 
JI=J-1 
J J  =  J-2 

W2=H23*CON*W( I , J ) 

W3=H120*CON*W< I , J ) 

C0W  =  W2* ( W(  1  »  J1  )-W(  I , JI  )  )-W3*( W (  l » J2  )-W(  I  ,  JJ  )  ) 

U  3  =  U ( I , J )*H1C 
U4=  U (  T  »  J )*H40 
V 3  =  V (  T  ,  J  )*H10 
V4=V<  I  ,  J  )*H40 

S A=  S (  T1,J)*(B1-U3)+S< I ,  J  1  )  *  <  B1 -V3 ) +S (  1 1  ,  J  )  *( B 1 +U3 ) 

1  +  S(  I , JI ) *(B1+V3) +S ( I  1 , J 1)  *<  B3+U4+V4)  +S(  I  I  ,  Jl  )  *<B3  +  V 

*  4  -  U  4  ) 

2  +S  (  I  I  ,  J  I  )*<  B3-U4-V4  )+S(  I  1  ,  J  1  )*(B3+U4-V4  )+S(I,JJ)*(B 

*  5-V4) 

3  +S(  IJ,  J)*<B5-U4)+S<  12, J )*( 85+U4) +S< I , J2 ) * ( 85  +  V4 ) -CO 

*  W 

SN=OMES*  (  S  A  -  S  (  I , J)  ) +S( I , J  ) 

IF( J.LT .NI .OR. J.GT.N)  GO  TO  14 
I F<  J  .EQ.NI  )  GO  TO  25 
S(  I  *N2)=-SN 
GO  TO  14 
25  SI  I »  N  3 )= - S  N 
14  SD  =  DABS( SN-SC  I  ,  J  )  ) 

0  I  =  DM AX  1 ( D I ,SD) 

HA=DABS ( SN ) 

REF=DMAX 1 ( REF , HA ) 

12  S (  I  ,  J)=SN 
SER=DI /REF 

WR  I T  E ( 6  »  8 1 )  SK »  SE  R »  OME  S 
IFISER.LE.  ERR.  OR.SK.GE. TIMS)  GO  TO  15 
S  K  =  S  K+ 1 
DO  27  J  =3 » (M 

S(2,J)=(S(3,J)*<B1-U3)+S(2, J1)*(B1-V3) 

1  +  S  (  2  »  J  I  )=MB1+V3)+S{3,J1 ) *( B3+U4+V4) +  S ( 3  »  J  I ) *{ B3 

*  +U4-V4) 

2  +S( 2, JJ )*( B5-V4)+S( 4, J )*( B5+U4)+S( 2, J2 ) *<  B5+V4) 

3  -COW  )/( 1D0+65-U4) 

27  S (M, J )  =  <  S(M, JI )*<  B1-V3 ) +S(M I , J )*( B1+U3) +S( M, JI)*( B1+V3 

*  ) 

1  +S(MI,J1)*{  B3  +  V4  -1)4  >  +  S  (  M  I  ,  JI  >  *  (  B3-U4-V4  )  +S  (  M  ,  J  J 

*  )*<B5-V4) 

2  +S(MJ, J)*(B5-U4)+S(M,J2)*(B5+V4) 

3  -COW  )/(  1DC  +  B5  +  U4) 

DO  2  8  I  =  3, MI 

COWI 2  =  CON*W ( I , 2 ) *H120*( W( I , 5) -6D0*W(  I ,4 )  +  18D0* 

1  W< I ,3)-10D0*W( 1,2)) 
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C0WI2=CQN*W{  1 , 2  )  *  (  H  8  0  *  W  (  I  ,  3  )  -H240*  W  {  1,4) 

1  -H144*W( I ,2 )+H1440*W( 1,5) ) 

28  S  (  I  »  2  )  =  (  S  (  I 1 ,2>*( B1-U3 )+S(  I , 3 ) *( 31 L-V3 )  +S  ( 1 1  *  2  )  *  ( B1+U3 ) 

1  +  S(  I  L v 3 ) *< B3+U4+V4)  +  S( I  I  * 3 ) * ( B3+V4-U4) +  S { I J , 2 ) * 

*  (B5-U4) 

2  +  S(  12,2  )*(B5+U4  )  +  S(  I  ,4)*{B5+V4) 

3  -COW  I  2 )/( 1D0+B5-V4) 

S<2,2)={ S(3,2)*{B1-U3)+S(2,3)*(B1-V3)+S{ 3, 3 ) * ( B3+U4+V4 

*  ) 

1  +S(4,2)*(B5+U4)+S(2»4)*{  B5  +  V4  ) 

2  -COW  I2)/(  l D0+2D0*B5-U4- V4 ) 

S ( M , 2 )  =  ( S ( M , 3 ) * ( B 1-V3 ) +  S ( M I , 2 ) * ( B 1 +U3 )  +  S ( MI , 3 ) *{ B3+ V4- 

*  U  4  ) 

1  +S( MJ,2 )*(B5-U4)+S(M,J)*(B5+V4) 

2  -COW  1 2  )  /  (  1 DO  +  2DO  *B5+U4- V4 ) 

GO  TO  16 

15  CONTINUE 

C  USING  STREAM  FUNCTION  TO  FIND  U  AND  V 

WRITE (6, 66) 

D  I  =  ODO 
REF  =0D0 
DO  17  J=2,N1 
J1=J+1 
J2= J+2 
JI=J-1 
J J= J-2 
DO  17  1=2, M 
11=1+1 
12=1+2 
11=1-1 
I  J=  1-2 

I F ( J  . EQ • 2 )  GU  TO  18 

UN=  (  S  (  I  ,  JJ  )-8D0*S<  I  ,  J  I  )+8D0*S<  I  ,  J  1  )  -S  (  I  ,  J  2  )  )  /  {  12DQ*H) 

GO  TO  19 

18  UN= ( S ( I , 5 ) -6 DO* S (  I,4)+18D0*S(  I  ,3)-10D0*S(  1,2)  )/(  12D0*H 

*  ) 

19  I F< I . EO. 2  )  GO  TO  20 
1 F { I .EQ.M)  GO  TO  21 

V  N  = ( S ( 12, J)-8D0*S(  II, J  )  +  8  D  0  *  S (  II,J)-S(IJ,J))/(  12D0*H) 

GO  TO  22 

20  VN= ( S ( 3  »  J ) -6D0*S ( 4, J )  +  1 8D0*S ( 3 , J ) - 10D0* S ( 2, J )  > / ( - 1200* 

*  H ) 

GO  TO  22 

21  VN= ( S ( MK, J ) -6D0*S( MJ ,J )  +  18D0*S { M I , J ) - 10 D0*$ ( M , J ) ) /  (  1 2D 

*  C*H ) 

22  VD=D ABS ( VN-V( I , J ) ) 

UD=DABS( UN-UI I , J ) ) 

DI =DM AX  1 ( DI ,UD, VD ) 

HA=D ABS ( UN ) 

HB=DABS ( VN ) 

REF  =  DM AX  1 ( REF , HA , HB ) 

U(  I,  J)=UN 
17  V  (  I  ,  J  )  =  V  N 
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VER=D I /REF 
WRITE(6,82)  K  ,  VE  R 

IF( VER.GE.VTEST.AND.VER.GE.5Q-2.AND.K.GE  .500 )  CALL  E  X  I 
*  T 

VTEST=VER 

ER1=ERR*2D0 

IF( VER.LE.ER1 .OR .K.GE.TIMU)  GO  TO  23 
K  =  K  + 1 
GO  TO  24 
23  CONTINUE 

THE  TEMPERATURE  PROFILE, T,(  THE  ENERGY  EQUATION  ) 

USING  THE  OBTAINED  W  AND  U»V 
T  K=  1 

WRITE (6, 67) 

11  D I  =  0  D  0 
R  EF  =  0D0 
DO  8  I  =  2  ,  M 
11=1+1 
11=1-1 
DO  8  J  =  2  »  N 1 
J1  =  J+  1 
J I = J-  1 

W1=H4*W<  I  ,  J) 

U2=U(  I , J )*H8*PR 
V2=V ( I , J )*H8*PR 

T  A  =  T {  I  1 , J) *< 25D-2-U2  ) +T ( I  I  , J ) * ( 25D-2+U2  ) 

I  +T ( I ,J1 )*( 25D-2-V2  )+T( I, JI )*(25D-2+V2  ) 

1  -  W  1 

TN=0MEE# ( T  A-l  II,J))+TII,J) 

I F ( J.LT.N.OR. J.GT.N)  GO  TO  9 
T  (  I  ,  N  2  )  =  T  N 

9  TD  =  DABS I  TM-TI  I  ,  J )  ) 

Di =DM  AX  1 ( D I , TO > 

H A=DABS ( TN ) 

R  EF  =  DMA X 1 ( REF , HA ) 

8  T  {  I  »  J  )  =  T  N 
TER=DI /REF 

WRITE (6, 80)  TKtTER»OMEE 
TF( TFR .LE.ERR.OR • TK.GE.TIME)  GO  TO  10 
TK=  TK+ 1 
GO  TO  11 
10  CONTINUE 

WRITE(6,49)  MfN,  ERR , GAMMA , PR , CON 

WRI T  E ( 6  »  50 ) 

WRITE  (5,51)  WK , WE  R , OMFM 

WRITE  (6,52)  ( ( W (  I, J), 1  =  2, Ml) , J  =  2 , N  1 ) 

WRITE  ( 6,53  ) 

WR I T  E ( 6  »  80 )  TK, TER, OWE 6 
WRITE(6,52)  (  (T(  T,J)  ,1=2, Ml)  ,J=2,N1) 

WR  I T  E  (  6 , 55  ) 

WRITE (6, 81)  SK  »  SE  R, OMES 

WRITE (6, 52) ( ( S(  I , J)  ,1=2  ,M1)  ,J  =  2,N1) 

WR I TE ( 6 , 54 ) 
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WRITE(6,8?)  K  ,  V  E  R 

WRITE<6,52) ( <U<  I, J)  ,I  =  2,M1)  ,J=2,N1) 

WRITE (6, 70) 

WRITE (6, 52) { <V( I, J) ,1=2, Ml) ,J=2,N1) 

C 

C 

C  THE  END  OF  ITERATIONS 

C 
C 
C 

c 

c 

C  CALCULATION  OF  FLOW  AND  HEAT  TRANSFER 

C 

c 

c 

C  WM  IS  MEAN  VELOCITY 

C  f M  IS  THE  MEAN  TEMP. 

C  BT  IS  THE  MIXED  MEAN  TEMP. 

C  TO  CALCULATE  MEAN,  MEAN  TEMP., MIXED  MEAN  TEMP.,  MEAN  S 

*  HEAR 

C  STRESS  AND  MEAN  TEMP.  GRADIENT  BY  APPLYING  SIMPSON  RUL 

*  E 

42  WM=0D0 
TM=CD0 
TW=ODO 

DO  34  1=2, M, 2 
WM=WM+2D0*W( I , N I ) 

TM=TM+2D0*T  (  I  ,  N 1  ) 

34  TW=TW+2DQ*W  (  I ,  N1  )  *T  (  I  ,  N 1  ) 

DO  35  1=3, MI, 2 

WM= WM+W (  I  , N 1 ) 

TM=TM+T ( I , Nl ) 

35  TW=TW+W(  I  ,N1 )*T (  I , Nl ) 

DO  36  J  =  2 ,  N  ,  2 

DO  37  1=2, M, 2 
WM= WM  +  W {  I , J  )*8D0 
TM=  TM+T (  I  , J ) *8D0 

37  TW=TW+W<  I  ,  J  )*T(  I  ,  J  )*8D0 
DO  36  I =3, MI, 2 

WM=  WM  +  W (  I ,  J  )  *4D0 
TM=TM+T (I , J ) *4D0 

36  TW=TW+W(  I  , J  )  *T(  I  , J  )  *4D0 
DO  38  J  =  3 » N I »  2 

DO  39  I  =  2 ,  M  ,  2 
WM=WM+W (  I  ,  J  )*4D0 
TM=TM+T (  I  ,  J )*4DC 
39  TW=TW+W( I , J)*T( I , J )*4D0 
DO  38  1=3, MI, 2 
WM= WM+W ( I , J ) *200 
TM  =  TM  +  T (  I  ,  J  )  *200 

38  TW=TW  +  W(  I  ,  J >*T(  I , J ) *2D0 
>jm  =  Ww*2D0/  (9D0*FM*FN) 
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TM=TM*?DO/( 9  00*1 M*FN ) 

T  W  =  T  W* 2 00 / ( 9D0* F  M*F  N ) 

WP I T  E ( 6*56) 

WRITE(6»68) 

00  500  J=1,N1 

SH  (  1  f  J)  =  (4DQ*W(  2  ,  J  )  -300* W(  3  ,  J  )  +4D0*W  (  4,  J  )  /  3D0 
1  -WC 5, J ) /4D0) /H 

TG( 1» J)=(40Q*T (2, J)-3D0*T(3, J ) +4D0*T ( 4 , J ) /3D0 
1  -T( 5,  J ) /4D0 ) /H 

SH( Ml  ,  J )=(4D0*W( M,  J  ) -3D0*W< MI , J ) +4D0*W( MJ , J ) /3D0 
1  -W( MK, J) /4D0) /H 

TG( Ml , J ) =( 4D0*T( M, J ) -300*1 <MI , J ) +4D0*T< MJ, J ) /300 
1  -T{ MK, J) /400 )/H 

500  CONTINUE 

DO  501  1=1, Ml 

SH(  I  ,  1)  =  {4D0*W(  I  ,2)-3D0*W(  I  ,  3)+4D0*W(  I  ,4)/3D0 
1  —  W I  1,5) /4D0 ) /H 

TG  <  I  ,  1  )  =  ( 4 0 0 * T (  I  , 2  )  -3D0*T (  I  ,3)+4D0*T< 1,4) /300 
1  -T (  1,5)  / 4 D 0 ) /H 

501  CONTINUE 

10  3  DO  104  J  = 1 , N 1 

WR  ITE(6, 57)  J  ,SH< 1 , J) ,SH (Ml , J)  ,TG( 1, J ) , TGCMl , J ) 

104  CONTINUE 
WRIT  E(6, 58) 

WR  I  T  E { 6 , 69 ) 

DO  105  1=1, Ml 

WRITE { 6, 59)  I , SH(  I  , 1  )  , TG< I ,  1) 

105  CONTINUE 

SHM=  SH (  1 ,  D+SHCMl,  1  ) +  SHC  1  ,N1  ) +SH  {  Ml  ,N1 ) 

T GM=TG(  1  ,  1  )+TG<  Ml ,  1  )+TG<  1  ,N1  )  +TG( Ml ,N 1 ) 

SUM=0D0 
SUT  =0D0 
DO  30  J  =  2  ,  N  ,  2 
SUM=SUM+SH{ 1  , J ) +SH ( Ml , J ) 

SU  T=  SUT  +  TG  (  1 »  J ) +TG( Ml , J ) 

30  CONTINUE 
SHM=SHM+4D0*SUM 
TGM=TGM+4D0*SUT 
SUM =0 DO 
SUT=0D0 

DO  31  J  =  3 , N  I  ,  2 
S  UM=  SUM  +  SH ( 1  , J ) +SH( Ml, J ) 

SUT  =  SUT  +  TG{ 1 , J ) +  TG{ Ml , J ) 

31  CONTINUE 

SHM=( SHM+2D0*SUM )  *H 
T  GM= ( TGM+2D0*SUT)*H 
SUM  =0  DO 
SUT  =  0D0 
DO  32  1=2, M, 2 
SUM=  SUM  +  SH ( 1,1) 

SUT=  SUT  +  TG ( 1,1) 

32  CONTINUE 
SHM=SHM+SUM*4D0*H 
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T  GM= T  GM  +$U 1 *4DQ*H 

SUM=ODO 

SUT=ODO 

00  33  1=3 , Ml, 2 
SUM=  SUM  +  SH (1,1) 

SUT=SUT+TG( 1,1) 

33  CONTINUE 

SHM=  (  SHM+2D0*H*SUM  )  *2D0*GAMM  A/  (  3  DO*  (  GAMMA-*- 1  DC  )**2  ) 
TGM=(TGM+2D0*H*SUT ) *2D0*GAMMA/ ( 300* ( GAMMA+ 1  DO ) **2 ) 

C  FRE  IS  THE  PRODUCT  OF  FRICTION  FACTOR  AND  REYIVOLDS  NO 

* 

C  NU  IS  THE  NUSSELT  NO. 

DK=C0N**5D-1*WM 

TB=PR*TW/WM 

FRE=?DO*SHM/WM 

FREM=?00/WM 

NU=WM*TGM/TW 

NUE=WM**2/ ( 4*  T  W ) 

WRI TE (6,60)  SHM, TGM, WM,TM,TW,TB 
WR  ITE ( 6,61  )  FRE, FRE M, NU,NUE ,DK 
CALL  PRINT  (WfT,S.»U,V) 

GO  TO  555 

45  FORMAT (3X, I  2, 2X, 13,4010.3,415) 

46  FORMAT! 5X, 025. 16, 2010. 3) 

49  FORMAT ( *  1  *  3 X ,2HM=, I2»3X,2HN=, 12,  5X  ,4HERR= ,010.3, 

1 5X , 6HGA  MMA  = ,015.4, 5X , 3HPR= , 0 10 • 3 , 5X , 4HC0N= , 0 1 0 • 3 ) 

50  FORMAT ( *0* 10X, 17HM0MENTUM  EQUATION) 

51  FORMAT ( 1  OX , 6HSWEE  P= , 13, 10X , 7HM A X ERR= , 0 1 3 . 6 , 10X,5HOMEM= 

*  ,013.6) 

52  FORMAT! IX,  8D13.5 ) 

53  FORMAT (• 0» 10X, 15HENERGY  EQUATION) 

54  FORMAT ( ' O'  1CX , 10HVELOC IT Y, U ) 

55  FORMAT (* 0* 10X, 15HSTREAM  FUNCTION) 

56  FORMAT ( • 0* 10X  *30HSHEARING  STRESS  OF  Y-R  IRECT ION  ,  12X , 

1 3 5H TEMPERATURE  GRADIENT  OF  Y-D  IR EOT  1  ON ) 

57  FORMAT ( 10X,2HJ=,I3,5X,,2(D13.5,7X),10X,2(D13.5,7X)) 

58  FORMAT (* 0* 1GX , 30HSHEAR I NG  STRESS  OF  X-R  I  REC T I ON , 1 2X , 

1 35HTEMPERATURE  GRADIENT  OF  X-DIRECTION) 

5  9  FORMAT! 10X, PHI =,  13, 5X, 01 3.5,7 X, 3 OX, 01 3. 5) 

60  FORMAT (  1  OX , 4HSHM = , 0 1 3 . 5 , 3X , 4HT GM= , 0 1 3 . 5 , 4X , 

l 3HWM= ,01 3. 5,4X, 3HTM= , 013.5 ,4X, 3HT W= ,D13. 5 ,4X, 3HT  B= , D13 

*  .5  ) 

61  FORMAT ( 10X , 4HFR E= , 013. 5 , 2X , 5HFREM= , 013. 5,4X , 3HNU=,D13. 

*  5 , 3  X  , 

14HNUE=,D13. 5,4X, 3HDK= ,013.5 ) 

64  FORMAT (  *  0* 10X,25HERRGR  OF  AXIAL  VELOCITY, W) 

66  FORMAT (  'O'  10X,29HERROR  OF  SECONDARY  FLOW,U  £  V) 

65  FORMAT (* O' 10X,24HERROR  OF  STREAM  FUNCTION) 

67  FORMAT (  *0* 10X, 19HERR0R  OF  TEMPER ATUR  ) 

68  FORMAT (20X,3HI=1 , 17X,4HI =M1 ,24X,3HI =1, 17X,4HI=M1 ) 

69  FORMAT<  20X , 3H J  =  1 ,45X,3HJ  =  1 ) 

70  FORMAT ( *0* I0X,10HVELOCITY,V) 

8  0  FORMAT ( 10X , 6HSWE  E  P= , I  3 , 1  OX , 7HM AX  ERR  =  , 0 1 3. 6 , 1 0 X , 5H0ME E  = 

*  ,013.6) 
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61  FORMAT ( 10X,6HSWEEP=* 13, 10X , 7HMAX ERR = , 01 3. 6, 10X,5HOMES 
*  ,  D  1  3 . 6  ) 

8  2  FORMAT (  1  OX, 6HS WEE P=, I  3 , 1 OX , 7HM AXERR= » 01 3.6 ) 

1000  STOP 
END 
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SUBROUTINE  PRINT  (M»T9StU«V) 

DOUBLE  PRECISION  W ( 8 1  * 43 ) » T ( 81 , 43 ) , S < 8 1 , 43) ,U( 81 ,43 ) »V 
*  (81,43) 

COMMON  M1,N3 

WRITE  (7,71)  ( (W(  I,J),I=1,M1)  »  J= 1 » N  3 ) 

WRITE(7 ,71 )  (  (T(  I  ,J)  ,1=1,  Ml)  ,  J  = 1 , N  3 ) 

WRITE(7»71) ( (S( I,J) ,1  =  1 ,Ml)  , J  =  1 , N3 ) 

WRITE(7,71)( (U(  I,J)  ,1  =  1, Ml)  ,J=1,N3) 

WR I T  E ( 7 , 7 1 ) (  (V(  I,J)  ,  I = 1 » M 1 )  ,J  =  1,N3) 

RETURN 

71  FORMAT ( 10  A8  ) 

end  I  1 
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1  MM+1 


Fig.  1  Coordinate  system  and  numerical  arid 
for  a  curved  rectangular  channel 
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Fig.  2  Effect  of  grid  size  on  flow  and  heat  transfer 
results  for  a  curved  square  channel  (y  =  1 )  with  Pr  =  0.73 


Computing  Time  (min.) 


Fig.  3  Effect  of  grid  size  on  computing  time  for  a  curved 
square  channel  (y  =  1 )  with  Pr  =  0.73 
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A  relation  between  error  e-,  and  computing  time 
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curved  square  channel  y  =  1  with  Pr  =  10 
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for  a  curved  square  channel  y  =  1  with  K  =  51.9  and  Pr  =  0.73 


140 

120 

100 

80 

60 

40 

20 

0 


Distribution  of  dimensionless  secondary  velocity 
ponents  (u,v)  in  two  directions  (X  and  Y)  passing 
hrough  the  center  of  circulation  for  a  curved 
square  channel  (y  =  1 )  with  K  as  a  parameter 
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Fig.  10  Dimensionless  axial  velocity  distribution  in  a  curved 
rectangular  channel  y  =  2  with  K  as  a  parameter 
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Fig.  11  Dimensionless  temperature  distribution  in  a  curved  rectangular 
channel  y  =  2  and  Pr  =  0.73  with  K  as  a  parameter 


0.50 


Fig.  12  Secondary  flow  streamlines  and  dimensionless  isothermals 
for  a  curved  rectangular  channel  y  =  2  with  K  =  58.8  and  Pr  =  0.73 
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g.  13  Distribution  of  dimensionless  secondary 
locity  components  (u,v)  in  two  directions  (X 
d  Y)  passing  through  the  center  of  circulation 
r  a  curved  rectangular  channel  (y  =  2)  with  K 
a  parameter 
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Fig.  14  Dimensionless  axial  velocity  distribution  in  a  curved 
rectangular  channel  y  =  5  with  K  as  a  parameter 
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Fig.  16  Secondary  flow  streamlines  and  dimensionless  isothermals 
for  a  curved  rectangular  channel  y  =  5  with  K  =  88.1  and  Pr  =  0.73 
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Fig.  17  Dimensionless  axial  velocity  distribution  in  a  curved 
rectangular  channel  y  =  0.5  with  K  as  a  parameter 
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Fig.  18  Dimensionless  temperature  distribution  in  a  curved  rectangular 
channel  y  =  0.5  and  Pr  =  0.73  with  K  as  a  parameter 


X/a 

Fig.  19  Secondary  flow  streamlines  and  dimensionless 
isothermals  for  a  curved  rectangular  channel 
y  -  0.5  with  K  =  103.4  and  Pr  -  0.73 
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Fig.  20  Distribution  of  dimensionless  secondary  velocity 
components  (u,v)  in  two  directions  (X  and  Y) 
passing  through  the  center  of  circulation  for  a 


curved  rectangular  channel  (y  =  0.5)  with  K  as  a  parameter 
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Fig.  21  Dimensionless  axial  velocity  distribution  in  a  curved 
rectangular  channel  y  =  0.2  with  K  as  a  parameter 
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Fig.  22  Dimensionless  temperature  distribution  in  a  curved  rectangular 
channel  y  =  0.2  and  Pr  =  0.73  with  K  as  a  parameter 
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Fig.  23  Secondary  flow  streamlines  and  dimensionless 
isothermals  for  a  curved  rectangular  channel 
y  =  0.2  with  K  =  79.3  and  Pr  =  0.73 


'  E\.0  *  -><l  bns  £.5E  =  >1  rt3N  S.O  *  r 


107 


■ 


:  ■ 


108 


»o 

o 


<o 

CM 


o 


>5 

X 


^01  x  6 


CD 

£ 

03 

cr 

<o 

XJ 

<u 

>  C\1 


£ 

O 

r— 

(J 

II 

ro 

£ 

£ 

Q_ 

•r— 

X3 

£ 

£ 

O 

fO 

•r— 

+J 

£ 

QJ 

-Q 

4-> 

•i— 

OJ 

£ 

E 

+-> 

03 

CO 

£ 

•r— 

03 

XJ 

Q. 

CL) 

03 

£ 

13 

CO 

+-> 

03 

US 

£ 

x: 

QJ 

Cl 

J£ 

E 

+■> 

<U 

•r— 

4-> 

5 

CO 

i — 

CO 

cu 

II 

T—“ 

£ 

o 

•1— 

r— 

CO 

QJ 

£ 

£ 

QJ 

£ 

E 

03 

■r— 

-£ 

Q 

O 

LD 

C\J 

• 

CO 

«r— 

Ll_ 

I 


109 


«o 

o 


«o  -Q 

04  V, 
O  >- 


o 


zO[*d 


t o 
o 


«o 

04 

o 


o 

o  ^ 
X 


»o 

04 

o 


•o 

o 

i 


a> 

s- 

tO 

cr 

to 

■a 

ai 

>  co 


s- 

o 

3 

1 — 

O 

II 

tO 

S- 

c 

Q_ 

■,— 

*o 

c 

c: 

o 

to 

•r— 

s- 

<D 

x> 

+J 

•i— 

<V 

s- 

E 

+-> 

to 

to 

s- 

•r— 

to 

~o 

CL 

Q) 

to 

s- 

3 

to 

4-> 

to 

03 

S- 

a> 

CL 

sz 

E 

+-> 

a> 

•p— 

+-> 

to 

, — 

to 

a> 

II 

r— 

c 

>- 

o 

•1 — 

to 

cu 

c 

c: 

d) 

c 

E 

to 

•i — 

JC 

a 

o 

to 

OJ 

• 

CD 

•r— 

Ll_ 

i 


no 


Ill 


001 


112 


113 


Fig.  30  9/ ( G ) q  versus  K  with  aspect  ratio  y  as  a  parameter  and  Pr  =  0.73 


Fig.  31  0/(0)q  versus  K  with  Prandtl  number  Pr  as  a  parameter  in  a  curved  square  channel 
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Table  1  Numerical  Results  for  a  Curved  Square  Channel  y 
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Table  3  Numerical  Results  for  Aspect  Ratio  y  =  0.2  with  Pr  =  0.73 
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